Low Complexity Approximate Bayesian Logistic Regression for Sparse Online Learning
Authors: Gil I. Shamir
Wojciech Szpankowski
Venue:  (2021)
Keywords: Bayesian methods
logistic regression
regret
online learning
MDL.

Abstract: Theoretical results show that Bayesian methods can achieve lower bounds on regret for online logistic regression.  In practice, however, such techniques may not be feasible especially for very large feature sets.  Various approximations that, for huge sparse feature sets, diminish the theoretical advantages, must be used.  Often, they apply stochastic gradient methods with hyper-parameters that must be tuned on some surrogate loss, defeating theoretical advantages of Bayesian methods.  The surrogate loss, defined to approximate the mixture, requires techniques as Monte Carlo sampling, increasing computations per example.  We propose low complexity analytical approximations for sparse online logistic and probit regressions.  Unlike variational inference and other methods, our methods use analytical closed forms, substantially lowering computations.  Unlike dense solutions, 
as Gaussian Mixtures, our methods allow for sparse problems with huge feature sets without increasing complexity.  With the analytical closed forms, there is also no need for applying stochastic gradient methods on surrogate losses, and for tuning and balancing learning and regularization hyper-parameters.  Empirical results top the performance of the more computationally involved methods.  Like such methods, our methods still reveal per feature and per example uncertainty measures.


Under review as a conference paper at ICLR 2021
Low COMPLEXITY APPROXIMATE BAYESIAN LOGIS- TIC REGRESSION FOR SPARSE ONLINE LEARNING
Anonymous authors Paper under double-blind review
ABSTRACT
Theoretical results show that Bayesian methods can achieve lower bounds on re- gret for online logistic regression. In practice, however, such techniques may not feasible especially for very large feature sets. Various approximations that, for huge sparse feature sets, diminish the theoretical advantages, must be used. Of- ten, they apply stochastic gradient methods with tuned some surrogate loss, defeating theoretical advantages of Bayesian meth- ods. The surrogate loss, defined to approximate the mixture, requires techniques Monte Carlo sampling, increasing computations per example. We propose low complexity analytical approximations for sparse online logistic and probit regres- sions. Unlike variational inference and other methods, our methods use analytical closed forms, substantially lowering computations. Unlike dense solutions, as Gaussian Mixtures, our methods allow for sparse problems with huge feature sets without increasing complexity. With the analytical closed forms, there is also no need for applying stochastic gradient methods on surrogate losses, and for tuning and balancing learning and regularization the performance of the more computationally involved methods. Like such meth- ods, our methods still reveal per feature and per example uncertainty measures.
that must be
Empirical results top
INTRODUCTION
We consider online (Bottou 1998 Chalev-Shwarde et al. 2011) binary logistic regression over series of rounds E 2I, T}. At round f, a sparse feature vector x; E [-,, with de < nonzero values, is revealed, and prediction for the label y, E {} must be generated. The dimension d can be huge (billions), but d, is usually tens or hundreds. Logistic regression used in huge portion of existing learning problems. It can be used to predict medical risk factors, to predict world phenomena, stock market movements, or click-through--ate in online advertising. The online sparse setup is also very common to these application areas, particularly, if predictions need to be streamed in real time as the model keeps updating from newly seen examples.
A prediction algorithm attempts maximize probabilities of the observed labels. Online meth- ods sequentially learn parameters for the d features. With stochastic gradient methods (Bottou 2010 Duchi et al. 2011), these are weights Wi,i associated with feature i E {1..... d} round f. Bayesian methods keep track of some distribution over the parameters, and assign an expected mixture probability to the generated prediction (Hoffman al. 2010 Opper & Winther [998). The overall objective is to maximize sequence likelihood probability, or to minimize its negative log- arithm. A benchmark measure of an algorithm's 's performance its regret, the excess loss it attains denoting transpose). A comparator w" that minimizes the cumulative loss can be picked to measure the regret relative to the best possible comparator in some space of parameter values.
over an algorithm that uses some fixed comparator values of w A (w1, ......, ..) (T
Kakade & Ng (205)); Foster et al. (2018): Shamir (2020) demonstrated that, in theory, Bayesian methods are capable to achieve regret, logarithmic with the horizon and linear with d, that even matches regret lower bounds for d = o(T'). Classical stochastic gradient methods are usually im- plemented as proper learning algorithms, that determine w, prior to observing X,, and are inferior in the worst-case (Hazan et al. [2014), although, in many cases depending the data, they still achieve logarithmic regret (Bach) 2010) 2014] [aach Moulines] 2013). Recent work (Jézéquel et al. 2020) demonstrated non- non-Bayesian improper gradient based algorithms with better regret.

---

Under review as a conference paper at ICLR 2021
Unfortunately, superiority of Bayesian methods diminishes by their intractabilityy A theoretically optimal prior has a diagonal covariance matrix, with each component either uniformly or Normally distributed with large variance. Effects of such prior cannot be maintained in practical online problems with a large sparse feature set, as the posterior of such a prior no longer has the properties subsequent prior. Gaussian approximations that rely on
of the prior, but must be maintained diagonalization of the covariance must be used. Neither normal nor the diagonal assumptions are true for the real posterior (even with diagonal prior). They thus lead to performance degradations. Diagonalization is similar to linearization in convex optimization, specifically for stochastic gradient descent (SGD) (Zinkevich 2000). It allows handling features independently, but limits performance.
Bayesian learning literature focused on applying such methods to predict posterior probabilities. and provide model (epistemic) uncertainty measurements (Bishop 2006, Dempster] 1968) Huelsenbeck & Ronquist 2001 Knill & Richards 1996). However, uncertainty of feature is, in fact, mostly a function of the number of examples in which it was present; a measure that can be tracked, not estimated. Methods, such as Variational Bayesian (VB) Inference (Bishop 2006 Blei et al. 2017 Drugowitsch| 2019; Drugowitsch et al. 2019; Ranganath et al.] 2014), track such measurements by matching the posterior. However, as demonstrated in Rissanen (1984) seminal work, minimizing gret is identical to uncertainty reduction, as regret is merely a different representation of uncertainty. Regret can be universally minimized over the possible parameter space through a good choice of prior. Hence, to minimize uncertainty, the role of an approximation is to preserve the effect of such a prior at least in the region of the distribution that dominates the ultimate posterior at the hori- zon T. This is a simpler problem than matching the posterior, and opens possibilities for simpler approximations that can lead results identical to those of heavy methods as VB.
VB methods are typically used offline to match tractable posterior to the true one by upper bound- ing overall loss. They are computationally involved, requiring either iterative techniques (Bishop 2006 Blei et al. [201]] Murphy 2012) like Expectation Maximization (EM) (Dempster et al.. 19777 Moon 1996); or Monte Carlo (MC) sampling, replacing analytical expectation by an empirical one over randomly drawn set. To converge, MC can be combined with gradient descent, either requir- ing heavy computations, or adapting stochastic gradient methods to update posteriors (Broderick) et al. 2013 Knowles 2015 Nguyen et al. 2017a b). For online problems, the posterior of an exam- ple is the prior of the subsequent one. To minimize losses, online VB must converge to the optimal approximation
every example. Otherwise, additional losses may be incurred, as the algorithm may not converge at each example, while the target posterior keeps moving with subsequent exam- ples. Moreover, combining methods that need tune hyper-parameters defeats the parameter free nature (Abdellaoui] 2000 Memahan & Streeter, 2012 Orabona & Pál 2016) of Bayesian methods.
Most Bayesian literature addressed the dense problem, where x, consists of mostly nonzero entries for every t, and the dimension d of the feature space is relatively small. Techniques, like Gaussian Mixtures (Herbrich al., 2003) Montuelle et al.) 2013; 2014; Rasmussen 2000; Reynolds et al. 2000 [ung] (0044,, that may use VB, usually apply matrix computations quadratic in d on the covariance matrix. In many practical problems, however, a very small feature subset present in each example. For categorical features, only one of the features in the vector is present at any example. Techniques, useful for the low dimensional dense problem, may thus not be practical.
We provide a simple analytical Bayesian method for online sparse logistic and probit regressions with closed form updates. We generalize the method also for dense multi- dimensional updates, if the problem is not completely sparse. Our results are first to study regret for Bayesian methods that are simple enough to be applied in practice. They provide an example to the connection between uncertainty and regret, and more broadly the Minimum Description Length (MDL) principle (Grunwald. 2007: Rissanen 1978a by 1984 1986, Shamir 2015) 2020). Empirical results demonstrate the advantages of our method over computationally involved methods and over other simpler approximations, both by achieving hetter regret and better loss on real data. As part of the algorithm, uncertainty measures are provided with no added complexity. We specifically demon- strate that it is sufficient to have an approximation focus on the location of the peak of the posterior and its curvature or value, which are most likely to dominate regret, instead of approximating the full posterior, which brings unnecessary complexity missing the real goal of preserving the effects of good prior. In fact, approximating the full posterior may eventually lead to poor generalization and overfitting by focusing too much on the tails of the posterior. Our approach directly approxi- mates the posterior, unlike VB methods that approximate by minimizing an upper bound on the loss. Finally, our approach leverages sparsity solve a sparse problem.
Paper Contributions:
2

---

Under review as a conference paper at ICLR 2021
Related Work: The simplest single dimensional online logistic regression problem (d 1 and F1. vt) was widely studied. Jefferys' prior, p(0) - 1/ (T is asymptotically
optimal (Clarke & Barron 1994, Xie & Barron 1997 2000 Drmota & Szpankowski, 2004)), It
can be expressed in terms log-odds weights w as p(w) - e" e"///// T (1 + eee)) Applying a mixture
leads the Krichevsky & Trofinov (1987) (KT) add-1/2 estimator Q (yile" = [n--1 (3)) +
0.5]/t, where (36) counts occurrences of ye- We use y¹ to express sequence from to Applying this prior in a Follow The Regularized Leader (FTRL) setting (McMahan 201 1
also
leads to the KT estimator. This raised the question whether regret optimality generalizes to large dimensions (McMahan & Streeter 2012). Hazan et al. (2014) showed that this was not the case for proper methods. The bounds, however, theoretically generalize for Bayesian methods (Kakade & Ng. 2005 Foster et al.) 2018) [Shamir 2020), with large variance Gaussian or uniform prior with diagonal covariance. Peaked priors fail, as for each feature in an example, other features provide a self excluding log-odds prior, that shifts the relation between the overall distribution and the feature weight. While wide priors are good theoretically, because of the intractability of the Bayesian mixture integrals, diagonal approximations that are used unfortunately degrade their effect.
Bayesian methods have been studied extensively for estimating posteriors and uncertainty (Bishop) 2006 Makowski et al.] 2002 [Sen & Stoffa 1996). There is ample literature researching such tech-
niques in deep networks (see, e.g., [Blundell et al. (2013); Hwang et al.] (2020) Kendall & Gal
(017) : Lakshminarayanan et al. (2017); Malinin & Gales (2018); Wilson (2020)) Most of the
work focuses on the ultimate posterior after the full training dataset has been visited. One attempts to leverage the uncertainty measurements aid in inference on unseen examples. Techniques like
expectation propagation (EP) (Minka 2001 2013) (see also (ishhop (006); Chu et al.] Cun- ningham et al.) (2011); Graepel et al. (2010)) and VB are used to generate estimates of the posterior. In a dense setup, where there is a relatively small number features (or units in a deep network), Gaussian Mixture models can also be learned, where a jointly Gaussian posterior is learned, usually with some kernel that is used to reduce the dimensionality of the parameters that are actually being trained. Such methods, however, do not fit the sparse online setup.
Variational methods are derived utilizing Jensen's 's inequality to upper bound loss of expectation by expectation of the negative logarithm of the product of the prior p(w) and data likelihood Ply" xJ_ w). Normalizing this joint by the expected label sequence probability gives the posterior P(w|x" xT yy).. Then, a posterior Q(-) with a desired form is matched by minimizing the KL diver- gence KT(Q)) P), which decomposes into expectation w.r.t. Q( over the loss on y² and K ((QlP))) between the approximated posterior and the true prior. The first term may require techniques like the iterative mean field approximation EM (Bishop 2006 Jaakkola & Jordan 1998), MC sampling to be approximated. Gradient methods can also minimize K LLQ] In the sparse setup, it is stan- dard to assume diagonal Q(-). In an online setting, the process can be iterated over the examples (or mini- mini-batchess, where the posterior at is the prior at I ++1 1. Computing the approximate posterior may be very expensive if done for every example. SGD can be used with MC sampling, but that would incur additional losses, as the posterior changes between successive examples. Like VB, HP minimizes the opposite divergence K (PIIQ) between the posterior and its approximate.
2 PRELIMINARIES
Let Pe(w) ) be the prior on the weights at round f, where we start by initializing some P1 (w). We will assume that p(-) is approximated by a diagonal covariance Gaussian, with means f4g.t and variances of, for component i at time t. Leveraging results in Kakade & Ng (2005); Foster al. (2018); Shamir (2020), if we restrict wi; [- B, B). uniform prior over this interval or a normal prior with standard deviation proportional to B can be picked. (To approximate a Dirichlet-1/2, O-mean normal prior with variance 2n can be used.) Observing sparse xx, the prediction for ye is given by
Pt P(v(x))) = / p(uelxe, w)pe(w)dw - / / w/x;
(1)
where for binary logistic regression. the probability of the label given the example and weights is given by the Sigmoid of the label weighted dot product of the example and weights
ply(x w)
A
e Sigma (yix7 w)
(2)
1 + exp (-Vix? w)
3

---

Under review as a conference paper at ICLR 2021
The expected prediction Pt in (1) marginalizes out the weights according to the prior pr(-) from the joint probability of w and yt- The prediction pt is a function also of all prior pairs sequence fx"=1, y "} through the prior pe(-). After observing yh, try to match (diagonal) posterior Q() to the weights that will equal the next round's prior
P+41(W) A Qe(w) ss p(w|x'; ,2) = p(ye Xx, w)px(w) p(ye Xr, w ())((w)
P())(x))
Pt
Using ST {xT= )" the logarithmic loss incurred by approximation Q(-) on the sequence of
predictions is L(Sr, Q) - Σi= log Let w" be some fixed comparator in the parameter
values' space. Then, the regret of approximation Q(-) relative to comparator w is given by
QUST, w*) A L(ST L(ST, w- Σ [logrt + log(1 + expl -xx w')]
(4)
The regret can measure the excess loss relative to the best possible w' comparator, if it is chosen.
3 MARGINALIZED BAYESIAN GAUSSIAN APPROXIMATION
In this section, we describe the proposed method. First, the Sigmoid is approximated by normal Cumulative Distribution Function (CDF). prediction for the label of the current example gen- erated shrinking the cumulative mean score function of the cumulative variance over all features. The main idea for updating feature distributions marginalizing away all other covariates for each feature in an example at given round, such that the mean and variance of the feature can be updated to match the location the peak and either its curvature or value for the true marginalized poste- rior. In Appendix B we demonstrate the same approach for Probil Regression. It follows the same steps, except that it does not require the initial approximation. Finally, Appendix D.I shows how similar approximation methodology can be used to apply simple multi-dimensional updates instead of marginalized one, the can be performed when sparsity limited.
Gaussian Approximation of Sigmoid: The relation between the logistic distribution and the Normal one was well studied in the statistics literature (see, e.g., [Bishop (2006); Murphy (2012). The Sigmoid function in (2) can be viewed CDF, which be approximated by normal CDF ((2) (The inverse of (-)) is the Probit function.) The derivative of the Sigmoid function is the O-mean Logistic Probability Density Function (PDF). Matching the PDFs, we have e" /(1 + ex x
1/V2To" exp (w /20²). This yields that the Sigmoid function can be approximated by a O-mean Gaussian CDF with variance 8/τ. Using the standard O-mean normal (-)) function, the argument
scaled by the inverse of the standard deviation / π //88 giving
Sigma(w)
(V V
(5)
Ite
More details about this approximation are in Appendix A
Approximation approach and some notation: With the diagonal and Gaussian assumptions, for each sparse example (with only de < d nonzero entries in xt), we can assume that we have a single normal random variable, whose mean is the weighted mean, and whose variance is the quadratically weighted sum of variances. Denote the example total weight, mean, and variance by
Wi
Pi.t:
o
(6)
(where the diagonalization assumption is important for the simplicity of the approximation of of).
Since we consider sparse problem, there is benefit to breaking the dependencies between features present in a given example and updating each independently. We can achieve that by marginaliz- ing the prior at t over all other features. Because we assume all features are jointly independent Gaussians, we can break the joint prior into a product of two components; one, the marginal of the feature, and the other the marginal of all other features together, i.e., the self excluding prior. To
4

---

Under review as a conference paper at ICLR 2021
match the posterior, then marginalize the latter, and match single dimensional posterior for each feature. We define the self excluding prior for feature at time t, its mean and variance as
d W, W
Σ jyh
p pa Ia.PE Lii - ox x oo (7)
and the single dimensional variable wt, we can compute in (T) replacing p(u, X,, w) in by normal CDF. Approximating this integral (see,
Prediction: With the probit approximation
e.g. Murphy (2012). Section 8.4.4.2, and Bishop 2006)) gives
yef4
Pr Sigma
(8)
+ 50
This result demonstrates how the prediction variance shrinks the prediction towards probability 0.5. Marginalization: Given the diagonalization assumption, the prior & be expressed as ps(w)
Ps,((w)) (w-)) where (-) is the prior on the self excluding prior of w.. Hence,
p(u,, w xx)) p(v: Xx,, (x:))
Marginalizing on w-; gives
p(u,, wx(x)) = (s,((
p(ye X (w--w pu.(w;))lw-dte
(10)
The inner integral, which marginalizes over w-; with its prior (w can be approximated by
(w (_
Iw-sd
exp
20 20-14
Sigma (2, + w- dw
2TO it
12
(w) - (_s,t)
IV Villznww + w dw
exp
2mT
20
2(2) [[V 8 (i + pe- + o- -1112) dz
[c) +
v U U1 Wi)
V1 + 802.1
Ve + wss)
Sigma
(11)
/1 No2..
Step (a) follows, again, from the approximation (5). For (b), we apply the change of variables v (w-: - (/-,,))//////) --_..//o-. where si-) the standard Gaussian PDF. The integral (b) gives with a = v ye (M + w;) and b² _ 0² to lead to (c). Finally, the approximation is used to go back from a Normal CDF to a Sigmoid in (d).
Vi+t2
Posterior: The posterior w; is given by plugging
into 10 normalizing Pe given in
1 Pt.
yelp i,t + 2e,
(w)) Qu.(ww.) = p(w, bo v')
(w,) Sigma
(12)
/ 00
The approximation on the right implies matching the current true posterior with the ith component of the approximate posterior Q(-). It can be simplified to
(w: (2,e+1)2
1
(_(wx
(wi m1,,))
(a in T.
exp
exp
Sigma
0i,+1
20,,-11
PioL.
20
1+ 002
(13)
5

---

Under review as a conference paper at ICLR 2021
Because the functional form of the posterior is not Gaussian, there multiple ways to fit Gaussian. We review alternatives in Appendix D However, we want to ensure that the regions of the true posterior much, as this will incur additional loss. is thus desirable to match the peak of the true posterior with the peak of the approximation. One method is match both the location and hight of the peak. The other, Laplace approximation (Bishop) 2006), matches the location and curvature at the peak. Both methods give the same approximate for pi,++1 but somewhat different one for T,4+1.
Approximations:
most likely to converge to at the horizon are not scaled down too
To give pi, we find w; that maximizes the r.h.s. of 13 or minimizes its negative logarithm. Let
(_u Us (u- i, + Fs.cHh,e11)
]
Yt (M-i, 1, + 2Lole)et1)
Pi,t+ - Sigma
1 ++ exp
exp
V V go
(14)
V+ go2re
except that Mi,e+1 replaces Piit and o²1,0 replaces of: Thus Paa+ is the prob- ability predicted for yt if we update May and shrink function of o. The minimization gives
be almost pa in 8
Pa.t1 Pa,t
Pi,t+).
(15)
1
Eq. (15) can be solved iteratively, where the Newton's method can be used, as described in Ap- pendix The solution for (1,,4+1 can also expressed in terms of the r generulized Lambert W
function Corlless al. 1996 Mezo & Bariez 2015).
Alternatielly, to avoid multiple iterations per update when using Newton's 's method, we can use a Taylor series approximation of Pi.c+ around 1 Pi,t, where
A Pt.
Sigma
Ut (u-i. + 2r.eM(e)))
Sigma
yeft V
(16)
+ 5 _
Like Pt,++ Pt, not Pt- Instead, it is the probability of ye as projected by the means of the weights f, shrunk as function of o²,,, instead of aff More importantly, it depens only on parameters before the update at + is applied, giving a closed form solution. Applying first order approximation we have
y f. (1 Pp.c)
Pi,11 Mi.t +
(17)
V 1 -o
- o2 (1-Pi,r) Pi.
If may be simpler store the precision 1/0 in which case, 17 may be easier compute by normalizing both numerator and denominator by of applying this normalization on the right term of the denominator. Second higher orders approximations necessary, as the first order one already gives identical performance to the iterative method.
also be applied, but may not be
After updating pi,+11, we can apply (13)
to update 0Tt-11 Plugging
15 in we solve for 0i,
Pi Pi,t
(14,441 m),)2²
P.
ye o
p((tt)] (18)
011+1
exp
20?
Pia++
exp
2(1 " "02
Alternatively (18). Laplace approximation can be used by finding the second derivative of the negative logarithm of the posterior, giving
V x.
of 2111
= 2
1+ oo
Pi.t (1 ps.a+)
(19)
The procedures described are summarized in Algorithm In Appendix D describe several different methods and approximations that can be used, including one that applies the same approxi- mation steps applied in this section without the marginalization. As empirical results, in the next section, show, however, there is no performance advantage to applying any of the more involved methods.

---

Under review as a conference paper at ICLR 2021
Algorithm Marginalized Bayesian Gaussian Approximation
1: procedure MARGINALIZED BAYESIAN GAUSSIAN APPROXIMATTONN(aramette:: Mo. og)
2:
Vi E E.,,.... ..... di Mi,1 + μo. of + To.
3: for (=1,2..... do
4:
Get
5:
Compute Ma. of with (6)
6:
Generate pt for y/t E (-1,) with 8).
7:
Observe Vt.
8:
for i: : * 0 do
9:
Compute Pi, Pi,c+ with [6 and [4), respectively, using Mi,c+1 Mi. for (4).
10:
Iterate on (15) and (14) with Newton method, or use (17) to update f1,++1
11:
Update of ofe+1 with either (8
12:
end for
end for
14:
end procedure
Mall Fastures Total Gearpie Logoddo
MAdS Fastures Total206, Per Gaurppie Wieen Lagodas
100
02
00
Oki 7500 L0KG
Gevn Neav M400600 60CO
Oa. 5100
203030
131000
6000
16:08
30200 400000
B0030
Figure 1: R2/ logt t vs. round f for different methods for randomly drawn d binary features, expected d./d /d features per example, and standard deviation of true log-odds noted in each graph.
NUMERICAL RESULTS
To benchmark regret of an algorithm, one needs the ground truth of real loss minimizing param- eters. On real data, the loss minimizing parameters are unknown. Furthermore, true benchmark datasets consist of non-stationary data, and the feature sets selected by model one trains may misspecify the "true" features of such a model. Thus real data may not give clean evaluation of benchmark dataset the end of this section, but
the proposed methods. We present results measure regret performance of different algorithms, we present results on synthetic data.
Synthetic Data:
We simulated data by setting d features with true log-odds weights that were drawn randomly with some prior (or with multiple priors, each governing subset of the features). At f, de & d features were selected in random from the of features, where different rules were used for different simulations to draw the de features. In the fully random case, we set a random fraction d parameter, and each feature was activated with probability a. We either used binary features Ti. (0. 1}, or also drew Fi, randomly in [0. 1]. For categorical features, the features were partitioned into categories, and for every example, one or preset number features from each category were randomly selected, with different types of randomness, including fast decaying long tail distributions. Let 0 the vector of true parameters, then, Pr(Y; = 1) was computed as Sigma(x? 0). The probability was used to randomly draw y/- We used Algorithm as well other algorithms, described Appendix D to sequentially predict yt and update. For gradient methods, used Stochastic Gradient Descent (SGD) with AdaGrad scheduling (Duchi al. 2011)
We ran grids of algorithm hyper- parameters for all algorithms to find optimal ones, and we show results for these optimal hyper-paramene for all algorithms. Since know the true weights,

---

Under review as a conference paper at ICLR 2021
Table 1: Runtime and regret coefficients rT = RT/ log for different algorithms on synthetic 200 features models with true log-odds std of 1 and algorithm parameters in Fig
Per Example Mean: 20
Per Example Mean: 40
1M Examples
10M Examples
IM Examples
10M Examples
Alg.
Time
FT
Time
Fr
Time
Fr
Time
rr
SGD
8s
117.6
58s
135.9
I0s
246.3
1:18m
219.88
ADF
6s
190.3
45s
540.3
7s
178.6
53s
388.51
DimGauss
9s
78.02
1:01m
82.54
14s
97.52
1:42m
109.16
Gauss
9s
77.66
1:05m
93.85
14s
91.5
1:32m
144.28
ApproxG
9s
77.66
1:05m
93.85
12s
91.5
1:32m
144.28
VB-100
3:47m
144.22
35:25m
649.57
7:21m
202.9
1:10:25h
1125.3
VB-1000
35:16m
85.07
1:11:01h
104.31
use them for comparator baseline w'. Curves show progressive validation regret. At round
measure the cumulative regret up given by R, - -EE=1 [log log(1 + exp(- w*)] Instead of showing R,, we plot R logt. If an algorithm has logarithmic regret, the normalized curve of the algorithm will converge to the constant. This methodology thus allows us to observe whether algorithm has logarithmic regret not. Results different variations (labeled by Gauss for updates with (15). and by Gauss Approx with (17) Labels also designate the prior used (o and po). and which variance approximation was used (G for (18)
shown for Algorithm I with its
and I. for (T9)). Reference results are shown for SGD, multi-dimensional Gaussian approximation update (DimGauss) described in Appendix D.I EP; using Assumed Density Filtering (ADF), and marginalized VB (VBApprox) described in Appendix D.3
Fig. D shows normalized regret for different true data configurations described the top of each graph with random binary features. More detailed results on multiple configurations are shown in Appendix E Unfortunately, unlike theoretical results (Shamir 2020). the prior has to fit the data for all methods of Gaussian approximation for good regret. This is true for any known practical Bayesian method, well as for the learning parameters of SGD. However, we choose a prior that matches the true prior, regret rates logarithmic in T, usually close to the lower bounds of 0.5 log T per parameter, are achieved on all experiments with Algorithm [ The same results are achieved whether we use the iterative version in (15) or its simpler approximation (17), and whether the variance is updated by (18) or (9). Unlike Algorithm | both, DimGauss and ADF appear to be optimized for priors that are different from the true one, and that depend on the fraction number of features that occur in each example. DimGauss requires larger prior with more features. As the sparsity is reduced, DimGauss with its optimal prior seems improve relative to Algorithm This expected, as the sparsity assumption becomes less valid. Algorithm| Toutperforms both ADF and SGD in all cases. We can find SGD each configuration, but are inferior to Algorithm 1 with increasing gaps with more active features.
that seem still exhibit logarithmic regret for
Table shows both execution runtimes and regret coefficients 'T Rr log T for the full simu- lations with the algorithms and configurations in Fig [ Simulations were machine and included synthetic data generation and similar outputs for all algorithms compared. they
single Ubuntu
DimGauss, VB, and ADF may have slight advantage, package, which is highly optimized for matrix operations. We show benchmarks for 10 and 107 examples. We observe rather equal runtimes for SGD, ADF, DimGauss and Algorithm[]] with slight advantage to ADF, which may be due to the highly optimized matrix operations. Good regret is obtained for both methods of Algorithm|] [ but DimGauss slowly improves over Algorithm I This is true because selecting 20 or 40 features out of 200 results in repetitions of co-occurrences as more data is observed, and DimGauss utilizes these co-occurrences. We show below and also observe from the results in Appendix D. however, that with higher degrees of sparsity, the DimGauss al- gorithm degrades, while Algorithm [ retains its advantages. Interestingly, the Newton method in the Gauss algorithm does not require runtime are matched between the two. This is attributed to the fact that the Newton method requires very few iterations converge. Boldfaced in the table are the poor runtime results of the VB meth-
implemented with the Eigen
time than its Taylor approximation, and both regret and
8

---

Under review as a conference paper at ICLR 2021
Gatgerkeat Tota 530660 Per Example Moare Logeddd S1D.2
Denchrunt Detract Prograshc valkatier Loglons DOC
Napoo PULLL
Onie Pd A.C.
4002
02
130
-14
Figure 2: Left: Re/ log I. vs. l. for synthetic data with millions of long tail features. Right: Loss relative to the best SGD on the Criteo benchmark dataset for multiple algorithms.
ods (shown with 100 and 1000 samples), whose worst-case complexity is Odd,NJT), where N is the number of samples and J is the maximal number of iterations. The table illustrates how runtime increases substantially relative to all other algorithms whose complexity is Od,TT). Regret only approaches that of Algorithm with 1000 samples, and is far inferior with 100 samples.
Fig. 2 (left) shows normalized regret of the different methods for a synthetic categorical model, with 22 categories. Category j.j 1,2,... consists of 2 features (a total of over 8M). For each example, single feature is randomly drawn from each category with Zipf distribution e/(n++ 1) 1 where c is normalizer and n is the feature index in its category. The true weights are randomly drawn as before, with standard deviation 2. This gives a long tail distribution over features, where for each category small indices are drawn more often, but many features from the long tail do oceur in examples. This models realistic sparse dataset. Algorithm [ outperforms other methods, and the best values of the DimGauss algorithm are inferior due to the sparsity.
The right graph in Fig. 2 giiees relative percent aggregate loss performance of the Bayesian algorithms relative to the best configuration found for an AdaGrad SGD on the Criteo datase[" We trained all algorithms on the over 45. examples in this dataset, which consists of 13 integer valued features, and 26 categorial features with different category counts. For each example, we generated prediction, computed its log loss on the label, and applied update. The aggregate log loss is a sum of data uncertainty and regret. The first term is equal for all algorithms and linear in the size of the data, where the second is the regret, which is sub-linear in the data size for good algorithm. Thus even small noticeable percent improvements imply possible substantial improvements of regret. Algorithm [ achieved with only linear models 0.465 progressive validation log loss, which is better, for example, than results reported also using deep networks in (Cheng et al. 2016). We observe advantages the Bayesian methods over SGD, where Algorithm [ was superior to all methods. The DimGauss method slowly degrades relative to the other methods due to the sparsity of some of the categorical features.
Criteo display advertising challenge benchmark Dataset:
CONCLUSIONS
We introduced a simple Bayesian mixture diagonal Gaussian approximation method based on marginalization for sparse online logistic regression and probit regression, that attempts to retain the affects of a good prior around the optimal values the weights. The method does not require the complexities of standard Bayesian methods, as VB, but was empirically shown to achieve regret rates as good and even better. With proper priors, empirical results were close to regret lower bounds, and superior to other Bayesian methods also measured with their best choices of priors. With the strong relation between regret and uncertainty, this approach gives good uncertainty estimates. The methodology was proposed for logistic regression, and extended for probit regression, but can be further extended to other settings. We also demonstrated a matching approach that performs high dimensional updates, and can be used for dense problems.
https: /ww. www.wwaaggag kaggle, com/ http////riiieee sp ay-ad- chal 1enge

---

Under review as a conference paper at ICLR 2021
REFERENCES
Mohammed Abdellaoui. Parameter-free elicitation of utility and probability weighting functions. Management science, 46(11):1497-1512, 2000.
Bach. Self-concordant analysis for logistic regression. Electronic Journal of Statistics, 4:384 1:384-414, 2010.
F. Bach. Adaptivity of averaged stochastic gradient descent to local strong convexity for logistic regression. Journal of Machine Leurning Reseurch, 15(1):595-627, 2014.
Bach and E. Moulines. con-strongly-conwee smooth stochastic approximation with convergence rate o(1/n). In Advances in neural information processing systems, pp- 773-781, 2013.
C. M. Bishop. Pattern Recognition and Machine Learning. Springer, New York, 2006.
David M Blei, Alp Kucukelbir, and Jon D McAuliffe. Variational inference: review for statisti- cians. Journal of the American statistical Association, 112(518):859-8777 2017.
Charles Blundell, Julien Cornebise, Koray Kavukcuoglu, and Daan Wierstra. Weight uncertainty in neural networks. arXiv preprint (0X.v:5505005224 2015.
L. Bottou. Large-scale machine learning with stochastic gradient descent. In COMPSTAT, 2010.
Léon Bottou. Online learning and stochastic approximations. On-line learning in neural networks, 17(9):142, 1998.
Tamara Broderick, Nicholas Boyd, Andre Wibisono, Ashia Wilson, and Michael I Jordan. Stream- ing variational bayes. In Advances in neural information processing systems, pp. 1727-1735, 2013.
Heng-Tze Cheng, Levent Koc, Jeremiah Harmsen, Tal Shaked, Tushar Chandra, Hrishi Aradhye, Glen Anderson, Greg Corrado, Wei Chai, Mustafa Ispir, et al. Wide & deep learning for recom- mender systems. Proceedings of the workshop on deep learning for recommender systems, pp. 7-10, 2016.
Wei Chu, Martin Zinkevich, Lihong Li, Achint Thomas, and Belle Tseng. Unbiased online active learning in data streams. In Proceedings of the 17th ACM SIGKDD international conference on Knowledge discovery and data mining, pp.
2011.
Bertrand S Clarke and Andrew R Barron. Jeffreys' prior asymptotically least favorable under entropy risk. Journal of Statistical planning and Inference, 41(1):37-60, 1994.
R. Corless, G. Gonnet, D. Hare, D. Jeffrey, and D. Knuth. On the lambert w function. Adv. Compu- 1996.
tational Mathematics,
John Cunningham, Philipp Hennig, and Simon Jacoste-Julien. Gaussian probabilities and expec- 2011.
tation propagation. arXiv preprint
Arthur P Dempster. A generalization of bayesian inference. Journal of the Royal Statistical Society:
Series (Methodological), 30(2):205-22 2322, 1968.
Arthur P Dempster, Nan M Laird, and Donald B Rubin. Maximum likelihood from incomplete data via the em algorithm. Journal of the Royal Statistical Society: Series (Methodological), 39(1): 1-22, 1977.
M. Drmota and W. Szpankowski. Precise minimax redundancy and regrets. IEEE Trans. Inf. Theory, IT-50:2686-2700 2004.
Drugowitsch. Variational bayesian inference for linear and
logistic
regression.
drXiv:131..543844.0 2019.
Jan Drugowitsch, André G Mendonça, Zachary F Mainen, and Alexandre Pouget. Learning optimal decisions with confidence. Proceedings of the National Academy of Sciences, 116(49):24872- 24880, 2019.
10

---

Under review as a conference paper at ICLR 2021
J. Duchi, E. Hazan, and Y. Singer. Adaptive subgradient methods for online learning and stochastic optimization. Journal of Machine Learning Research, 12:2121-2159, Feb. 2011.
Dylan J Foster, Satyen Kale, Haipeng Luo andMehryar Mohri, and Karthik Sridharan. Logistic regression: The importance of being improper. In COLT Conference on Learning Theory, 2018.
Thore Graepel, Joaquin Quinonero Candela, Thomas Borchert, and Ralf Herbrich. Web-scale bayesian click-through rate prediction for sponsored search advertising in microsoft's bing search engine. Omnipress, 2010.
P. D. Grunwald. The Minimum Description Length Principle. MIT Press, 2007.
E. Hazan, T. Koren, and K. Y. Levy. Logistic regression: Tight bounds for stochastic and online optimization. In The 27th Conference on Learning Theory, COLT 2014, pp. 197-209. MIT press, 2014.
Ralf Herbrich, Neil D Lawrence, and Matthias Seeger. Fast sparse gaussian process methods: The informative vector machine. In Advances in neural information processing systems, pp. 625-632, 2003.
Matthew Hoffman, Francis R. Bach, and David M. Blei. Online learning for latent dirichlet alloca- tion. In Advances in Neural Information Processing Systems 23, pp. 856-864. Curran Associates, Inc., 2010.
John P Huelsenbeck and Fredrik Ronquist. Mrbayes: Bayesian inference of phylogenetic trees. Bioinformatics, 17(8):754-755, 2001.
Seong Jae Hwang, Ronak Mehta, Hyunwoo Kim, Sterling Johnson, and Vikas Singh. Sampling free uncertainty estimation in gated recurrent units with applications to nonmative mod- cling in neuroimaging. In Uncertainty in Artificial Intelligence, pp. 809-819 819. PMLR, 2020.
Tommi S Jaakkola and Michael I Jordan. Improving the mean field approximation via the use of mixture distributions. In Learning in graphical models, pp. 163-173. Springer, 1998.
Rémi Jézéquel, Pierre Gaillard, and Alessandro Rudi. Efficient improper learning for online logistic regression. In COLT Conference on Learning Theory, volume 125 of Proceedings of Machine Learning Research, pp. 2085-2108. PMLR, 09-12 Jul 2020.
M. Kakade and A. Y. Ng. Online bounds for bayesian algorithms. In NIPS, pp. 641-648, 2005.
Alex Kendall and Yarin Gal. What uncertainties do we need in bayesian deep learning for computer vision? Advances in neural information processing systems, pp. 5574-5584, 2017.
David Knill and Whitman Richards. Perception as Bayesian inference. Cambridge University Press, 1996.
David Knowles. Stochastic gradient variational bayes for gamma approximating distributions. arXiv preprint arXiv:1509.0000 2015.
R. E. Krichevsky and V. K. Trofimov. he performance of universal encoding. IEEE Trans. Inform. Theory, 11-2772):19972077, Mar. 1981.
Balaji Lakshminarayanan, Alexander Pritzel, and Charles Blundell. Simple and scalable predic- tive uncertainty estimation using deep ensembles. Advances in neural information processing systems, pp. 6402-6413, 2017.
David Makowski, Daniel Wallach, and Marie Tremblay. Using bayesian approach parameter estimation; comparison of the glue and meme methods. Agronomie, 22(2):191-203. 2002.
Andrey Malinin and Mark Gales. Predictive uncertainty estimation via prior networks. In Advances in Neural Information Processing Systems, pp. 7047-7058, 2018.
Brendan Mcmahan and Matthew Streeter. No-regret algorithms for unconstrained online convex optimization. In Advances in neural information processing systems, pp. 2402-2410, 2012.

---

Under review as a conference paper at ICLR 2021
H. B. McMahan. hebbowthh-rggullriiee--eeeee and mirror descent: Equivalence theorems and L1 regularization. In AISTATS, 2011.
H. B. McMahan and M. J. Streeter. Open problem: Better bounds for online logistic regression. In Journal of Machine Learning Seeearch-Procedding Track. 23, 2012.
1. Mezo and A. Baricz. On the generalization of the lambert w function. / 2015.
Thomas P Minka. Expectation propagation for approximate bayesian inference. arXiv preprint arXiv:1301.2294 2013.
Thomas Peter Minka. A family of algorithms for approximate Bayesian inference. PhD thesis, Massachusetts Institute of Technology, 2001.
Lucie Montuelle, Erwan Le Pennec, and Serge Cohen. Gaussian mixture regression model with logistic weights, penalized maximum likelihood approach. arXiv preprint 1r4iv:1304.2666 2013.
Lucie Montuelle, Erwan Le Pennec, et al. Mixture of gaussian regressions model with logistic weights, a penalized maximum likelihood approach. Electronic Journal of Statistics, 8(1):1661- 1695, 2014.
Todd K Moon. The expectation-maxxiizztiie algorithm. IEEE Signal processing magazine, 13(6): 47-60, 1996.
K. P. Murphy. Machine Learning A Probabilistic Perspective. MIT Press, Cambridge, Mas- sachusetts, 2012.
Cuong V Nguyen, Thang D Bui, Yingzhen Li, and Richard E Turner. Online variational bayesian inference: Algorithms for sparse gaussian processes and theoretical bounds. 2017a
Cuong V Nguyen, Yingzhen Li, Thang D Bui, and Richard Turner. Variational continual learning. arXiv preprint arXiv: 1710: 2017b.
Manfred Opper and Ole Winther. A bayesian approach to on-line learning. On-line learning in neural networks, pp. 363-378, 1998.
Francesco Orabona and Dávid Pál. Coin betting and parameter-free online learning. Advances in Neural Information Processing Systems, pp. 577-585, 2016.
Rajesh Ranganath, Sean Gerrish, and David Blei. Black Box Variational Inference. volume 33 of Proceedings of Machine Learning Research, pp. 814-822. PMLR, 2014.
Carl Edward Rasmussen. The infinite gaussian mixture model. In Advances in neural information processing systems, pp. 554-560, 2000.
Douglas A Reynolds, Thomas Quatieri, and Robert B Dunn. Speaker verification using adapted gaussian mixture models. Digital signal processing. 10(1-3):19-41 2000.
Rissanen. Modeling by the shortest data description. Automatica, 14:465-44 471, 1978a.
J. Rissanen. Minimax codes for finite alphabets. EEE Trans. Inform. Theory, IT-24(3):389-392 1978b.
J. Rissanen. Universal coding, information, prediction, and estimation. IEEE Trans. Inform. Theory, (T-30(4):629-636, Jul. 1984.
Rissanen. Stochastic complexity and modeling. The Annals of Statistics, 14:1080-1100, 1986.
Mrinal K Sen and Paul L Stoffa. Bayesian inference, gibbs' sampler and uncertainty estimation in geophysical inversion 1. Geophysical Prospecting, 44(2):3133350, 1996.
Shai Shalev-Shward Shalee---hwrzeeeec et al. Online learning and online convex optimization. Foundations and trends in Machine Learning, 1(7):107-19, 201
12

---

Under review as a conference paper at ICLR 2021
Gil I. Shamir. Minimum description length (MDL) regularization for online learning. In Proceedings of the 1st International Workshop on Feature Extraction: Modern Questions and Challenges at NIPS 2015, pp. 260-276, 2015.
Gil I Shamir. Logistic regression regret: What's the catch? In COLT Conference on Learning Theory, volume 125 of Proceedings of Machine Learning Research, pp. 3296-3319, 2020.
Jack Sherman and Winifred Morrison. Adjustment of an inverse matrix corresponding to a change in one element of given matrix. The Annals of Mathematical Statistics, 21(1):124-127, 1950.
Hsi Guang Sung, Gaussian mixture regression and classification. PhD thesis, 2004.
Andrew Gordon Wilson. The case for bayesian deep learning. arXiv preprint arXiv:2001.10995. 2020.
Q. Xie and A. Barron. Minimax redundancy for the class of memoryless sources. IEEE Trans. Information Theory. pp. 647-657, 1997.
Q Xie and A. Barron. Asymptotic minimax regret for data compression, gambling, and prediction. IEEE Trans. Information Theory, 46:431-445, 2000.
M. Zinkevich. Online convex programming and generalized infinitesimal gradient ascent. In ICML, pp. 928 928-9228868 936, 2003.
A RELATION BETWEEN GAUSSIAN AND SIGMOID
The Sigmoid function, which converts log-odds to probability is very close in shape to the Gaussian Cumulative Distribution Function (CDF) (22,, as well established in the statistics litcrature (see, e.g.. Bishop (2066) Murphy (2012)). The derivative of the Sigmoid function is given by
dSigma(w) dw
(20)
(1 e
and equals the PDF of O-mean Logistic distribution. We can approximate the logistic PDF by Gaussian by matching the PDFs,
w²
exp
(21)
+ em
V2TO
20
Matching the distributions at w = yields o = V8/T,
1
T
w²
+
exp
16
v 2m
exp
2 5
(22)
Thus, can approximate the Sigmoid with 0-mean Gaussian CDF with variance 8/T, giving (5)
It remains to demonstrate that the PDFs (and CDFs) are close to each other not only at the peak. [eeuueBdemonstreee the approximation of the Logistic distribution (left) and the Sigmoid function (middle) by a normal PDF and a normal CDF both with variance 8/T,, respectively. The green curve shows the differences between the logistic/Sigmoid and the normal, which are also plotted on the right plot at larger scale. The magnitude of the difference between the Sigmoid and the normal CDF is bounded by 0.02 over the whole region. The differences appear asymmetric around the origin, and are substantially small at 1/3 standard deviation from the origin, or less. While they can accumulate over multiple examples, it appears that the most probable scenario that positive and negative differences over multiple examples cancel each other. Furthermore, the motivation of a Bayesian method is to converge toward a peaked point mass, at the loss minimizing value of the parameter. As the variance is narrowed closer to such point mass, the approximation tends to exist in the flat region around the origin, where the difference between the Sigmoid and the normal CDF very small.
13

---

Under review as a conference paper at ICLR 2021
Figure 3: Left: Logistic and Normal (00, 8/m² distributions and their differences. Middle: Sigmoid and Normal CDF (0, 8//)) and their differences. Right: Differences between logistic and normal (PDFs and CDFs).
B PROBIT REGRESSION
In this appendix, we show the derivation of the method proposed in this paper for Probit Regression, where, in a similar manner to 2 the predicted label probability with weight vector w, label yt. and covariates is given by the normal CDF
yix! 1 2
a= (-2
xxx7 [ a(a)da 2 (yex? w)
p(uelxe; w)
exp
do
(23)
where, recall, o(-) and ((-) are the standard Gaussian (normal) PDF and CDF, respectively. While for logistic regression, we used a Gaussian approximation to obtain analytical expressions for the prediction in (8) and the marginalization integral in (TT), for probit regression, these no longer approximations. For the posterior, we will still apply a Gaussian and a diagonal approximations, as in the derivations based on (13).
The approach for probit regression is similar to the one described in Section 3 for logistic regression. For each feature track the mean Mi,e and the variance of, for the ith feature.
Prediction:
For example t, we use (6) to compute the total weight w,, its mean μe and variance of Eq. 7 gives the self excluding weights, their means, and their variances. Similarly to (8), using the approximate normal prior at t, we can derive the label prediction for y/:
W/ mx)
p P(v((x))
exp
+ (yix? w)
dwv
2TO?
20
(w M1)
rewe
exp
exp
dz dw
/ /2mo²
207
22T
2
2
Ve (e(orot)))
1
2a
exp
exp
dz
dv
2m
2
v 2T
2
Φ(v) T [velore + 14] de
yift o?
+
(24)
For (a), we use the definition of w, in 6 Step (b) follows from substituting v = (w, p+)/on
Step (c) identified the integrands as a product of the standardized A 1) normal PDF multiplied
by standardized normal CDF (a,v + me). This integral gives normal CDF &b VI+ for
y/A, and 8 o. of leading to (d).
Marginalization: Following the marginalization steps in Section 3 we can express the joint proba- bility of weight w, and label yı conditioned on the covariates x, and marginalized over all the other
14

---

Under review as a conference paper at ICLR 2021
nonzero covariates example in (10 by
yr ye + k(
p(yr, ws/x)) ps,r(u;))
(25)
1 o=..
where we use the steps of (IT), excluding the approximations, derive (5
Posterior: The posterior on wy is given as in 12p normalizing p(ye, wy Xe) by Ps from (24
Ye ye(u-i, i,t + vi,w))
Pi,+1((u;)) v p(ux 3²) = Pt
Pr.((w))
(26)
1 -ii
This posterior can now be matched by a normal posterior Qi,((w:) as 13
Approximation: Here, we follow the Laplace approximation applied in Section 3 All other meth- ods mentioned the paper are also possible. To find (,,,e-1 that minimizes negative logarithm of the r.h.s. of (26). define, similarly to (14
Ve(p-ir + Fiefiie11)
Zes+
(27)
11+021
as the probit score, which serves as the argument the normal CDF, where the ith mean has been updated, but all other means have not. We can now express the update of the ith mean by
MiP.,107 ((2,,)))
p4,+1
Pi,
11+021
((2,,)))
(28)
Eq. (28) is a similar update for probit regression to that of (15) for logistic regression, where the ratio (i, 1D (2i,4+) replaces 1 Pi,++ (and the scaling of the self excluding variance is unnecessary).
As 15 (28) must be solved iteratively because the term of ((zz,++)/ (zz,z)+++)) is function of (4,,411 through the definition of 2,,+- As in Section 3 we can use a first order Taylor approximation of ((2,,c)) 3/b((()++))) around its value for Mi. Similarly to (16). we define
WelU-i, + Ii,((ir))
Zi
(29)
V V1+021.
which is the score before update of all means Mi, but unlike the one used compute Pt, normal- ized by the ith self excluding variance o2 instead of o?. With some algebra, this gives single operation update, similar to that in (17), given by
pi Pig
(30)
1 (1+ 1 v o
o(er.) 4(x..)
25,t +
d(a,.)) (c...))
The term 2,, + ((zz,t)(((zz)) in the denominator replaces Pr, in the logistic regression update equation.
Taking the second derivative of the negative logarithm of the posterior and approximating 1/0 by it, gives a single operation update of the variance, similarly to 19
-
{ oop
((E,))))
((EL++)
- ((2,14))
21,t
$ ((
(31)
C MEAN AND VARIANCE UPDATES
Eq. [5) gives an update for the mean M44+1 that cannot be solved in closed form. This is beacuse Piac+ a function of pi.+1. However, the update is easily solvable with few iterations of Newton's
15

---

Under review as a conference paper at ICLR 2021
method.
We start by plugging Mit++ Mest for iteration E We follow by computing pi+t+ with The solution for 1i,+ is the value of wy that minimizes the negative logarithm of the r.h.s.
(f=0)
(0)
14).
of 12 where (r,((w;) is Gaussian with mean Pi, and variance o..-- At iteration we can compute
with (14. using (. The gradient of the negative logarithm of the posterior w.r.t. wy at
Pi,++
W/ = Mi,+1 is given by
(f-1)
(8-1) Mi,t+1
Mi,t
(1 (1-Pi pi,t+ ((-11)
9t
(32)
o,t
11 + g
The second derivative w.r.t. w; is given by
(e-1)
(1- pl,tt ((=11)
(1-p((=))
h
- i,t
+
1 50
Pist
(33)
Then, Mi is updated by
(6) (f-1) gr Mi,t+1
M1,t1 h
(34)
The process terminates when the difference |i,++ Mite1 is smaller than some threshold, specified number of iterations had already executed. Note that the second derivative (33) (the
(E-1)
Newton Hessian step) gives the same expression as the update the precision (inverse variance) in
(9)).
Eq. (88) gives update to the variance. It is interesting to interpret some observations from this update. For Ii, 2 0 Pr.++ 2 This leads to ye - pa,e) > 0 because for every Vt, Mi,c+1 has move away from ui, with the sign of yt- This implies that we add positive term from the log-odds converted to Pt to those converted into Pi,s+ and also apply less shrinkage using o² instead of of yielding this claim. Hence, the coefficient of og outside the exponential term in (18) is upper bounded by 1. The argument of the exponential term is always nonnegative, which implies that the exponential term lower bounded by 1. However, the self-excluding variance o²,, large, it makes the argument small. Similarly, if Pi,++ is large relative to pt it makes the whole expression smaller. These two observations imply that the self excluding prior is less certain, or if the self excluding prior has different belief about the label from the one observed overall, the uncertainty of the current feature reduces more, because it is deemed responsible for the observation Vr- On the other hand, if the opposite holds, i.e., either the uncertainty of the self excluding prior is low, or the self excluding prior agrees more with the observed label yt, then, feature matters less for the observation, and therefore, its uncertainty is not reduced much.
D OTHER METHODS
In this appendix, we describe updates we perform with other Gaussian approximation methods. We can update the mean vector and covariance of de features that occur at f, without marginalization, discarding off-diagonal covariance terms. This approach, that otherwise uses similar approximations to the ones we used in Section is described in Section D.I where some of the approximation steps are novel this paper.
Alternatively, one can choose any two points on the true posterior and match Q(-) on these. Least squares can be used with several points that represent a region of w;, for which the posterior at 7 likely to have most mass. EP and VB can be applied, discussed in (occioon[DD2]nndl]
respectively.
For the latter, specifically, single dimensional VB on the marginalized posterior can be applied to
minimize K KL (w:)) pp,, [x", )))) where p(w: xx v²) is the true posterior on the r.h.s. of (2) given by p(y,, wo (x(xx)/pt where p(gt: wi (x)xx) is in (10). This still requires expensive Monte Carlo estimation of expectations with an iterative Newton method. Without marginalization, however, similar VB approach would require even more expensive Monte Carlo sampling or iterative mean field approximation EM converge all components.
16

---

Under review as a conference paper at ICLR 2021
D.1 MULTI-DIMENSIONA GAUSSIAN APPROXIMATION
Instead of marginalizing on all other features to update wi for which Fi,t * O, we can apply multi- dimensional update on all features for which Fi, A round I. Such updates will enhance cor- relation between these features, and may be a better fit to problems in which such correlation is expected. For this update, we assume that the true posterior consists of product between prior with a diagonal covariance matrix and a Sigmoid, and we apply Lapace approximation to obtain new mean vector and covariance. With some abuse of notation, let all values t consist only of the de nonzero components of x, Let Eg be the diagonal covariance matrix, with diagonal elements of !.- Let u, be the estimated mean vector at t. Then, the true posterior at is given by
exp 2 (w u, T (w u,)
1
p( x
(35)
Pt.
(2 de |E
+ exp yyxx w|
Similarly to
14). define
Pc+ _ Sigma (VIX? u
(36)
as the probability of y computed with weights after they had been updated (and this time with no shrinkage), where u is the updated vector of means. Then, with Laplace approximation, taking the value of the mean vector that maximizes the posterior, the mean can be updated as in (13) by
u u, + E. (1 ppt)).
(37)
This is, again, an equation that must be solved either numerically, or using methods such as Newton's method. Again, we can assign u ue and apply (36) on it to obtain p:_.. Then, at iteration E,
gr = Er -
(u (6-1)
uz) - yext (1 -1-ptt"
P++
(f-1))
(38)
and
H(0) = 2t + yixxx T. (l-1) (1 ((-1)).
Pit
(39)
Then, 1441 is updated by
(t-1) Hir) ge.
(40)
Termination is either when the update on all components of 1 is less than some threshold, or after a set number of iterations. Inverting the Hessian H also gives the updated covariance Σ++1, whose diagonal elements can be now used to update of,41, if we apply the algorithm for sparse problem, where it is infcasible to store all covariances.
Instead of updating Hig, we can keep track of its inverse HiT). and there is no need to invert the covariance matrix E,, With the diagonal form of E,, all operations can implemented with linear complexity in de using the Sherman & Morrison (950) formula, which simplifies matrix inversions for special matrices. For our specific need here, if is some matrix, a some constant, and x some vector, then, the Sherman-Morrison formula is
--xx" A
(A + oxxx) = A-1 1 oxTA-Ix
(41)
Substituting A = D, x = xx, and a = -yiph
(8-1) "1 (1-ptt 1),
p2+ p(--1)),,
we update (0) inverting 39
As in the marginalization method described in Section 3 we can avoid the iterative Newton method
with a first order Taylor approximation of - Pc+ around - Pr. where p. is defined in an analogy to (16) as
p. Sigma (vxx7 u.)
(42)
as the un-shrunk prediction of yt at round (which is different from pe, which is shrunk by the variance). The approximation leads to the following set of equations to update both u,+1 and E,+1. For simplification, define
V. Σ,x,
(43)
Then, temporarily update Dx, using Sherman- Morrison formula, to
- (1 T
vi - v.)
Σ = E
= Σ,
(44)
1 1 5 Σxxe
1 I+ vř Fx)xx Ve
17

---

Under review as a conference paper at ICLR 2021
Since Σ, is diagonal, the transpose the last term of the numerator in the first equality is unneces- sary. The second equality gives vector multiplications. showing that the complexity is linear in the dimension of the vectors dz (This is true also to the computation of v, when Σ, is diagonal.) Next, u can be updated
Ur+1 = + ve -.:)- P.)Ex-IX.
(45)
Now, we can update P++ in (6). using U++1. and use it to update Ee+1 using Sherman-Morison,
vipe+ (1
E E E,
(46)
+ vi Pt+ pe+ (+++)x' x Vi.
In the sparse case, we can now take the terms of the diagonal of E++1 to update oi of the nonzero covariates at round
Finally, it may be simpler to update the precision matrix He+1 E instead of the covariance Ee+1. Specifically, if multiple updates are performed in mini batch, the update applied to the covariance cannot be applied additively. However, additive updates on the precision are valid. Thus the updates in (44) and (46) can be replaced by
H = H. + v.Pc -x))xxx
(47)
and
HA+1 = H. + v. Pc+ -1- xx
(48)
respectively. To update U141, we still need to invert H4+1. We can use (44) if update was applied to single round only. If mini-batch update additively applied multiple updates at once in (47) the updated Hi+1 must be inverted to obtain E,+1-
The multi-dimensional approach described in this section can be applied to sparse problems, but also to dense problems. In the dense case, the operation in (43) is no longer linear in de, as the covariance matrix not necessarily diagonal. The use of Sherman-Morrin Morrison formula, however, to invert the precision and covariance, still applies and lowers the complexity of the approach. In the sparse problem, however, this approach may try to force correlations that are not there, that are then ignored. As empirical results suggest, it may not be as good as the marginalization approach because of that. Furthermore, unlike the marginalizatton approach in Section 3 which achieves best performance if the true prior matches the used to initialize the algorithm, empirical results demonstrate that the best performances are obtained with priors that are diffcrent from the true one with the multi-dimensione method when applied sparse problems.
D.2 EXPECTATION PROPAGATION ASSUMED DENSITY FILTERING
Instead of minimizing the divergence between the approximate Q and the true posterior, can the expectation propagation approach, as proposed in Minaa (00T), which essentially minimizes the opposite KL divergence, and attempts match the first two moments. More details can be found in Minka 2001).
D.3 MARGINALIZED VARIATIONAL BAYES
Instead of using Laplace approximation or matching the location of the peak of the true posterior and the estimated one together with either its curvature its value, we can apply full VB, by matching the approximate posterior Q with the true one through minimizing the KL divergence K L[QI) be- tween Q, the approximate posterior, and the true posterior. This requires either iterative approaches, such as mean field approximation EMM. or Monte Carlo sampling in order to approximate expectation over yet unknown Q. One can apply this approach in d, dimensions as the Laplace approximation in Subsection D.I] However, due to the inability to separate the covariates (in the Sigmoid), we would require a power set of samples. use 4d samples. This can be infeasible and complex if there are a large number d, of nonzero covariates. Instead, we can use VB only in the ith dimension for each feature separately together with the marginalization proposed in Section 3 This can be done by matching the approximate posterior (w,,(w;) with the posterior p(w:(x* xx', we obtained on the r.h.s. (12 on feature after marginalized on all other features.
we use N samples per dimension, this approach would
18

---

Under review as a conference paper at ICLR 2021
The KL divergence can be decomposed into three terms; the KL divergence between Q and the prior Ps.t, the contribution of conditioning the posterior on the probability Pt predicted for yt and the log loss (negative log likelihood) term, emerging from the Sigmoid.
[ (4- + Di. Wi)
=
KL(Qi,llPA.) + EQ.. log + EQ
log 1 exp
V o
EQ.. -EQ..
ouPP log
(Wi-Her -H.+11+11
(We mlle)
0i,+1
207,411
20
Y1 (M-1, + 2..)Wo)
log
1
1 exp
+ 3022.0
log
ot
1
1
0itt1 + 1h,t1 + ph,,- 2pr
01t+1
2
20
H
Y/ (u-i,t + Fz.e W.)
Equa
log
exp
+ exp
(49)
V1 V 31
All expectations are w.r.t. Qi,-- The KI. term be computed in closed form, giving the second and third equalities.
The last term the r.h.s. of (49 cannot be analytically computed without knowledge of the posterior Qi at t. Instead, we use Monte Carlo, by drawing N samples Sj (O 1), letting Wi,
Mie+1 SSoi,t4 With known 14,4 and 0i,411 we now approximate the expectation term in 49 as
N log
ye (u-2a ia The (M6,e1 + 8801.++1))
Σ 3=1
exp
where sj the jthh randomly drawn sample. Unfortunately, fi,+1 and o must updated in this step, and are not known. This requires, again, an iterative update using Newton's 's method. Similarly (14 in Section 3 we need to define a prediction Piji++ for which the prior (time t) means and variances used for all covariates except the ith one, and the updated 14,44 and of used for the ith mean and variance, respectively. This time, however, this prediction is defined N times, uniquely for each of the jth samples
Y1 (M-ie 24. (14,44 + s30+
Pijac+ Sigma
(50)
V + _
Then, the updated mean satisfies
1 N
Mi,
Σ (1 Pi.a+)) N
(51)
I 50
However, order to compute both Mi,++1 and 0 we need to apply Newton's method, optimizing (0)
(0) Mi.++ and 02 together. We start with Mi,++1 = f4,t and 0 di.t- For simplicity, we omit
the iteration number (e) from the notation. The following should be read as updates at iteration &
(8-1) which use Pi.j.t for the updates. For simplicity, let yrTa,e/ / / (* 18)8 Then, the joint
A
gradient w.r.t. Pi. and oe, is given by
21 N j=1
(1 Pi...+)
N
g
ΣN 8y Pi.j.t)
(52)
o
19

---

Under review as a conference paper at ICLR 2021
The joint Hessian is given by
oi I + N Eラ= Pi.je- Pi.jet)
oi
Ei 8jPi -P(g.4))
N
H
de E== spPi.jtt (1 - Pi.jet)
oun + of + a n Σソ= sfpt,j.4(1 Pi.jet))
N
(53)
Finally, at iteration l
M Mi++1
- -
(54)
(2)
1,0+1
1,1+1
As before, the update terminates if the differences between the values of two iterations are less than some threshold, or after set number of iterations.
As observed, this method requires J) operations for a single update (and ((d,NT) oper- ations overall), where the set number of Newton iterations. Empirical results demonstrate that even with N as large as 1000, results are not as good as those with the marginalization with Laplace approximation, presented in Section We note that one can use two dimensional first order Taylor approximation on Pi.j.e+ around Ps.j.t. which is defined similarly to Pt.j.t+, with the exception of using He,t and Fg, instead of ps,e+1 and 0i 1. obtain an approximation for updating pi,++1 1 and 01t++1 as in (17). The approximation should be made for every j. It will, however, not require the O(J) operations of the Newton method. The complexity still of O(N) factor greater than that of the method in Section 3
E ADDITIONAL EMPIRICAL RESULTS
In this appendix, we bring collection of simulation results demonstrating the performance of Algo- rithm and the other methods in different settings. all the simulations we performed we observed that Algorithm]] [ with proper prior for the setting, consistently gives regret close to the lower bound (e.g., 0.5 0.5 T cost for each unknown parameter). The other methods, while in some cases exhibit performance close that of Algorithm I
fail to do consistently in all conditions. The marginalized VB approach requires very large complexity as shown in Table and Fig. approach the performance of Algorithm
Fig. gives a two dimensional grid of varying d, and varying true feature weights. In all these simulations, Algorithm [ with prior matched the true one, gives minimal regret curves. While SGD seems, with the right choices of parameters, to approach logarithmic regret, its regret is larger, and increases with the true parameters' variance and number of active features de. Multi-dimensiooal Gaussian approximation appears to be mis- mis-calibrated on the prior. and depending on the feature density de tends to become better only with much larger priors. EP ADF gives larger regrets and while appearing reasonable with higher feature density, seems to be inferior with lower feature densities. Similar results are obtained when replacing the data generation model with models with uniform priors on the feature weights with various ranges.
Fig. 5 shows curves for simulations with different models. On the left, features are nonbinary, and on the right, the model consists of order of magnitude more features, where in average order of magnitude more features occur in each example. In both cases, Algorithm [ persists with similar regret rates, whereas other algorithms exhibit larger regrets. With large number of features, both the multi dimensional Gaussian approximation and SGD have larger regrets, although with even lower prior/ / learning rates slightly better regret may be possible.
Fig. demmonstrreess curves of the various algorithms for categorical features, where in each example fixed of features are selected from each category. On the right, the selection gives exponentially decaying distribution over the features, so that some features are selected very often while others rarely occur. Again, regret rates behave similarly for A [n and SGD. The multi-dimmesional Gaussian approximation completely breaks in this setting. This is because it hypothesizes correla- tions in its updates with features that rarely reoccur.
20

---

Under review as a conference paper at ICLR 2021
Figure 4: Normalized R,/ log vs. round for various methods with randomly drawn binary fea- tures, with d, expected de/d, and standard deviation of true log-odds noted in each graph. Graphs shown for d = 200, Eld:] E (5, 20, 40} and true log-odds std in 21,
Tesal
STD
203
PLL
Npraa
102
170
PAR
Figure 5: Normalized Re/ log vs. round & for different algorithms and different data generation models. On the left, model with d = 200 features, out of which in average d, 40 occur in an example, and weight standard deviation 2, with nonbinary feature values uniform in [0, 1]. Right: models with d 2000, average d, 200, data generation standard deviation of 2, with binary features.
21

---

Under review as a conference paper at ICLR 2021
Categoricat Tale 2900 Per Crappl Muar: 3, Legodda STD: 2
Calegericll Tolal0000, Per Grample Maan: Logodds STD:
Darda PCAD
000
110
60003 I000 1ed0
20000 40000 6X004
8000
Figure 6: Normalized by R,/ logt t round t for categorical models and different algorithms. Left: 2000 features 10 categories 200 features each, where in each example 5 features from each category are present, with true weight standard deviation of 2. On the right: A similar setting, except that in each category features are selected with a long tail exponential distribution (prioritizing a few features and rarely selecting others).
22


Reviews:
[{"review_id": "gPkhpm7ZHkQ", "reviewer": "AnonReviewer1", "summary": "", "questions": "", "limitations": "", "rating": 4, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "PAPER SUMMARY\n\nThis paper proposes a new posterior approximation scheme for probabilistic logistic regression. In the online learning setup, the proposed approximation leads to a closed-form parameter update. In particular, the key ideas here include: \n\n(1) A probit approximation of the sigmoid likelihood in Eq. (5). This, along with the assumed statistical independence among individual features and another integral approximation trick used in (Murphy, 2012), leads to a posterior representation in forms of a normalized product between a sigmoid and a normal PDF. \n\n(2) Matching the location and density of its mode to that of another Gaussian surrogate to yield the closed-form update that expresses the parameter at round (t + 1) in terms of its counterpart at round (t) analytically in Eq. (15)\n\nNOVELTY & SIGNIFICANCE\n\nIf I understand correctly, the main arguments in favor of the proposed method are: (a) the iterative update derived from Eq. (15) only involves a single sample at round (t) which is more efficient than VB's or EM's stochastic gradient updates that require re-iterating over multiple mini-batches of past data; and (b) the new approximation scheme performs empirically better than existing method, and matches state-of-the-art regret lower-bound.\n\nGiven the current form of the manuscript, these arguments, however, are not very convincing to me as I elaborated in the followings:\n\nFirst, regarding (a), while it is true that VB's or EM's offline updates are more expensive than the update derived from Eq. (15), it seems to me they can still be made more efficient in online setting with very minor modifications. For instance, in both VB and EM, the lower-bound function often factorizes additively across data point so if we cache all gradient computation in step (t), the computation of new gradient in step (t + 1) also involves only the latest data sample which is likely not more expensive than iteratively solving Eq. (15) -- I would like to hear the authors' thoughts on this.  \n\nSecond, the claim that the new approximation scheme performs better than existing methods is not well-supported. In this regard, I wonder why the improved loss is only demonstrated in a very controlled, synthetic experiment settings. One could always argue that as far as only synthetic experiments are involved, there might exist other settings for which we might see the opposite. To me, positive results on a single set of synthetic experiments often do not speak conclusively to whether one method would be more useful than another in practice -- a better way to demonstrate this is to also evaluate the proposed method on real-world benchmark datasets (especially those extracted from practical domains that were mentioned in the introduction)     \n\nFurthermore, the evaluation does not seem to include all the most relevant baselines (more on this later) as well as processing time comparison. Given that one of the key contribution claims here is the computation advantage over EM's and VB's, both theoretical complexity and empirical demonstration of averaged running time should be provided. \n\nIn addition, I also find part of the positioning of this paper somewhat misleading, especially when it criticizes stochastic gradient methods for incurring additional losses while the target posterior keeps moving with subsequent samples. Isn't this also true for the proposed method? While Eq. (15) is analytic, solving it for the parameter update is clearly not analytically tractable which also requires iterative update; and then even the resulting approximated posterior has to be matched again to another Gaussian surrogate. Lastly, the conclusion surprisingly summarizes that with proper prior, the method \"matches regret lower-bound\" -- this seems like an overclaim as there is no theoretical analysis in the paper to back this up.\n\nTECHNICAL SOUNDNESS\n\nI have gone through the derivations and I do not spot any serious issues, but there are parts that seem to lack justification. For instance, could the authors please detail the steps that lead to Eq. (5)? Also, it will be even better if the authors could provide plots to visually demonstrate how the derivative of sigmoid is close to the scaled normal CDF.\n\nI also find it a bit strange that this paper highlights the sparse aspect of data as an issue but in the end, it appears more to be a necessary setting that justifies the assumed diagonal form of the variance matrix, which is a necessity for most of the technical derivations in this paper.  \n\nCLARITY\n\nWhile the paper is sufficiently clear, I find the introduction is a bit too long. It could have made its points using much less space. In the current form, it kind of distracts the readers from those key points by throwing in many subsidiary positionings that are better discussed & organized in related work. \n\nThe claim on analytic update is also a bit exaggerating as solving Eq. (15) for an update equation is after all not tractable and we will have to resort to iterative method such as Newton's or Taylor series expansion anyway. \n\nEXPERIMENT\n\nOn the experimental evaluation, my key concern (as mentioned above) is the lack of evaluation of real-world dataset, and in addition, some of the more recent methods such as (Nguyen, 2017a; 2017b) were not included in the baseline.\n\nThe presentation of the evaluation is also problematic when the results seem to be collected on only 1 single run. As a standard practice, please re-run the same experiment multiple times and report the margin of error. It is also unclear why the Y-axis is t/log(t) -- first of all, do you mean r_t /log(t) instead?; and secondly, why do we divide it by log(t) but not t? Aren't we ultimately interested in the limit of the average regret r_t/t?\n\nFurthermore, despite the claim that this method is proposed for the setting with huge feature sets (with billions of dimensions), the synthetic experiment is mostly around 200; and up to 2000, which is a bit disappointing. On the same note, despite having a main claim of computational efficiency, there is not a single experiment that showcase this advantage over existing methods.\n\nREVIEW SUMMARY  \n\nThe paper aims to develop a posterior approximation scheme for sparse online logistic regression. While the key idea is to achieve faster computation time and better regret rate than existing methods, this point has not been demonstrated sufficiently due to a complete lack of complexity analysis and experiments on large-scale, real-world data with more recent benchmark. The paper also appears to overclaim on the analytic tractability of its update as I have pointed out above. \n\nNote: This is only my preliminary assessment and there is of course the possibility that I might miss some important points. In that case, I am looking forward to receiving detailed clarifications from the authors.", "extra_scores": {}, "extra_text": {}}, {"review_id": "ZWikB0AbKdl", "reviewer": "AnonReviewer4", "summary": "", "questions": "", "limitations": "", "rating": 4, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "The authors propose a low complexity approximation method with closed analytic forms for doing logistic regression in the sparse, online setting. They first introduce the marginalized bayesian gaussian approximation approach, which essentially replaces the sigmoidal with a gaussian. They then give approximate expressions for the prediction and marginalization terms, as well as ways to approximate/update the posterior. \n\nThis is an interesting approach. The pros are clear: closed-form analytic updates that attempt to retain the core advantages of a Bayesian approach. \n\nThe authors should address the following questions/concerns:  \t \n\t \t\t\nAssumptions: The diagonal Gaussian assumption appears to be quite strong. Are there comments/extensions on how to incorporate this for correlated/collinear covariates? \n\nApproximation: The approximation is based on a Gaussian approximation to the sigmoidal. Are there any theoretical/quantitative guarantees on the approximation error? Right now the approximation is just a heuristic.  \n\nExperiments: since one of the main proposed advantage of this scheme is that it is computationally faster, it would be nice to have experiments/charts showing how much faster (in actual experiments) does this method run when compared to competitors like VB or a Gibbs sampler. Also, there is a Newton’s method step in the algorithm. Could the authors comment on how that affects run time? \n\nMotivation: If you are using a Gaussian to approximate the Sigmoidal, since you are using a Gaussian anyways, why not just use an online probit regression then, instead of the logistic? \n\nOverall, I think this is a potentially useful proposal that is interesting, but I think a) more experiments needs to be done with competitors to illustrate not just regret, but runtime comparisons b) there needs to be more theoretical rigor as to how an approximation performs and what the approximation error/bounds/guarantees are. \n\nI also looked at the citations more closely, and I discovered that there are some mis-citations: \nFor example, you cited Terry Anderson. The theory and practice of online learning. Athabasca University Press, 2008, which upon further inspection is about online (as in over-the-internet) learning, not online learning in the machine learning sense. \n", "extra_scores": {}, "extra_text": {}}, {"review_id": "78Z1BNcz4zA", "reviewer": "AnonReviewer3", "summary": "", "questions": "", "limitations": "", "rating": 4, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "This paper proposes an algorithm  for online logistic regression  based on analytical formulas that approximate the Bayesian predictive posterior. These approximations are based on assuming a diagonal covariance Gaussian form for the posterior at each iteration that is optimized to fit the true posterior. Two alternatives are proposed for the optimization: one based on Newton's method and the other based on a Taylor series approximation, which experimentally yields similar results. Thanks to these closed formula approximations, the resulting algorithm has a constant low cost per observation and  is particularly suitable for sparse high-dimensional scenarios. It is empirically shown that the algorithm achieves a known regret lower bound on synthetic data when the true prior on the weights used to generate the data is provided to the algorithm.\n\nStrong points:\n- The introduction properly positions the work\n- A low-complexity algorithm for online prediction suitable for sparse high-dimensional features is proposed\n\nWeak points:\n- No theoretical guarantees on the regret\n- Regret is strongly dependent on the prior\n- No experiments on real datasets\n\n \nThis paper proposes a low-complexity approximated bayesian algorithm that aims at attaining known regret lower bounds that have been recently shown [Shamir 2020] to be achievable by exact Bayesian methods. Unfortunately, no theoretical guarantees are provided on the impact of such approximations and the empirical performance is reported to be strongly affected by the choice of the prior. In particular, it is reported that the prior has to fit the data for good regret. This strong dependence on the prior makes it of limited practical interest in real life online scenarios. Therefore, I think that the paper is not strong enough and recommend its rejection.\n\nDetailed comments:\n\nSection 3 could be better structured and a \"road map\" at the beginning of the section would help the reader to quickly understand the different steps involved in the approximation. \n\nIf I understood correctly, in order to build a synthetic dataset, \n- you sample the weights $w^\\star$ using some prior (a gaussian with mean 0 and STD specified in the title of each plot).\\\n- At time $t$, you sample random features $x_{i,t}$ using different distributions, depending on the type.\\\n- You add up these features weighted by $w^\\star$ to obtain the true log-odds that is used to define the probability distribution (through the sigmoid) from which $y_t$ is sampled.\n\nThen, the following sentence is confusing:\n\"The “true” log-odds of the features drawn at t were added up, weighted by $x_{i,t}$, and converted to probability with the Sigmoid function.\"\n\nIn the following sentence you use $\\tau$ that doesn't appear in eq. 4:\n\"Results describe a time series on (4), where summation is on $\\tau$ up to example $t$ at time $t$.\"\n\n\"We demonstrate that it is sufficient to approximate the component of the posterior which will dominate\nat the horizon, matching it by a diagonal Gaussian approximation ,\"\nHow is this demonstrated?\n\nEq 11: (d) should be with $\\approx$ instead of $=$\n\nSince nonparametric Bayesian methods are mentioned in the abstract, it would be interesting to see a comparison against online nonparametric Bayesian methods.\n\nIn my opinion, the title should say \"online prediction\" instead of \"online learning\", since the latter does not necessarily imply predicting at each round and cumulating the loss. \n\nTypos:\\\n -- mostly a function *of* the number\\\n -- Jeffery's prior -> Jeffreys prior\n\n\n=====POST-REBUTTAL COMMENTS======== \n\nI thank the authors for the response and the efforts in the updated draft, in particular with respect to the new experiment on the Criteo dataset. \nRegarding the following arguments:\n\n1- *\"The dependence of regret on the prior applies to any Bayesian method, in the same way that tuning the learning rate for gradient methods is necessary.\"*\n \nI agree that, for batch methods, the prior can be tuned to the data. Nevertheless, in an online prediction setting, the prior needs to be chosen before observing the data and once you start predicting you can't go back and change the prior for the already predicted sequence. \n\n2- *\"To the best of our knowledge, none of the literature that proposed practical (online) Bayesian methods proposed regret bounds or even studied regret for these methods\"*\n\nRegarding both arguments, there is recent literature that proposes a practical nonparametric online prediction method with regret guarantees for a fixed given prior:\n\n*Lhéritier, Alix, and Frederic Cazals. \"Low-Complexity Nonparametric Bayesian Online Prediction with Universal Guarantees.\" Advances in Neural Information Processing Systems. 2019.*\n\nTherefore, unfortunately, my first two concerns (theoretical guarantees on the regret and strong dependence on the prior) remain so I retain my original decision.", "extra_scores": {}, "extra_text": {}}, {"review_id": "5p--WCT6f6O", "reviewer": "AnonReviewer2", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "##########################################################################\nSummary:\n\nThe paper proposes an algorithm for learning the parameters of a logistic regression model in an online setting. The proposed algorithm is based on two approximations: the posterior at iteration t over the model parameters is assumed to be multivariate Gaussian distribution with a diagonal covariance matrix and the logistic/sigmoid function is approximated by the CDF of the normal distribution. The numerical results show the usefulness of the proposed algorithm.\n\n\n##########################################################################\nReason for score:\n\nMy overall impression is that the proposed algorithm is of sufficient interest to warrant the acceptance of the paper. \n\n\n##########################################################################\nPros:\n\n1. The paper proposes analytical expressions for updating the parameters of the logistic regression model in an online setting. The computation cost for updating the parameter is reduced.\n\n2. The paper is overall well written and the authors describe in detail each step in the development of the proposed algorithm. \n\n##########################################################################\nCons:\n\n1. For me it is not clear why you settled on a logistic regression problem instead of probit regression and present the logistic regression as an extension. From my understanding, the keystones for the proposed algorithm are the approximation of the prior as a multivariate Gaussian with diagonal covariance and the marginalization step (which itself is due to the Gaussian assumption). From my point of view, the introduction of the extra approximation layer due to the use of the logistic regression is detrimental to the easy understanding of the paper. I believe it would have been better to just introduce a section that describes the modifications to make in order to apply the proposed algorithm to the logistic regression problem. \n\n2. In section 1 on page 2 you make a statement about the role of an approximation in minimizing the the uncertainty. For me it is not clear how the proposed algorithm addresses this statement. Could you please elaborate on the subject? Also, on page 3 you say that it is sufficient to match the component of the posterior that dominates at the horizon. What do you exactly mean by the component of the posterior? \n\n\n##########################################################################\nMiscellaneous:\n\n1. It's welcoming that the proposed algorithm has very good results, however I'm a bit surprised that the results are that good. My wariness comes from the fact that in equation 13 you match an approximation to the posterior distribution, which itself was obtained using an approximation. I would expect the errors due to the approximations to accumulate as more and more data items are processed. However, the regret seems to remain steady as you increase the number of items seen. Do you have an explanation as to why the errors due to the approximations do not accumulate?\n\n2. There are some typos here and there in the article, however I believe you can easily detect them with the help of a spell checker.", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: This paper proposes a new approximate algorithm for Bayesian logistic regression in the online setting. The primary approximation involved in the algorithm is the use of a diagonal Gaussian approximation. (A probably more minor approximation is approximating the sigmoid with a Gaussian.) The main discussion focused on two issues: Firstly, there was some sentiment that the paper lacked theoretical guarantees. Second, there were concerns about the experimental results. I feel that it is not a serious flaw that the paper lacks a theoretical regret bound. Given the current state of algorithms for this problem practical algorithms remain very much of interest. However, the general sentiment of reviewers was that the experimental results were not as strong (or as numerous) as would be hoped. For an algorithm without a theoretical regret bound, I do agree that stronger empirical evidence would be expected. This was partially addressed in a revision but still I agree with the consensus that more extensive numerical evidence should be expected, and for that reason I am recommending rejection.

Finally, I'll mention some other issues that I view as not counting substantially against the paper. Firstly is the dependence on the prior. Here I am in agreement with the authors that this is an aspect shared by all Bayesian methods. This issue is a (valid) argument about the value of all Bayesian methods, but not one I think we will resolve here. Second, there were suggestions from the reviewers about improvements that could be made to the baseline methods. Here I don't feel that it's fair that we ask the authors to make novel improvements to other algorithms, unless those improvements are very "obvious".
Decision: Reject