Gradient-based training of Gaussian Mixture Models for High-Dimensional Streaming Data
Authors: Alexander Gepperth
Benedikt Pfülb
Venue:  (2021)
Keywords: Gaussian Mixture Models
Stochastic Gradient Descent
Unsupervised Representation Learning
Continual Learning

Abstract: We present an approach for efficiently training Gaussian Mixture Models by SGD on non-stationary, high-dimensional streaming data.
Our training scheme does not require data-driven parameter initialization (e.g., k-means) and has the ability to process high-dimensional samples without numerical problems.
Furthermore, the approach allows mini-batch sizes as low as 1, typical for streaming-data settings, and it is possible to react and adapt to changes in data statistics (concept drift/shift) without catastrophic forgetting.
Major problems in such streaming-data settings are undesirable local optima during early training phases and numerical instabilities due to high data dimensionalities.%, and catastrophic forgetting when encountering concept drift.
We introduce an adaptive annealing procedure to address the first problem,%, which additionally plays a decisive role in controlling the \acp{GMM}' reaction to concept drift.
whereas numerical instabilities are eliminated by using an exponential-free approximation to the standard \ac{GMM} log-likelihood.
Experiments on a variety of visual and non-visual benchmarks show that our SGD approach can be trained completely without, for instance, k-means based centroid initialization, and compares favorably to sEM, an online variant of EM.

Under review as a conference paper at ICLR 2021
GRADIENT-BASED TRAINING OF GAUSSIAN MIXTURE MODELS FOR HIGH-DIMENSIONAL STREAMING DATA
Anonymous authors Paper under double-blind review
ABSTRACT
We present an approach for efficiently training Gaussian Mixture Models (GMMs) by Stochastic Gradient Descent (SGD) with non-stationary, high-dimensional streaming data. Our training scheme does not require data-driven parameter ini- tialization (e.g., k-means) and has the ability to process high-dimensional samples without numerical problems. Furthermore, the approach allows mini-bater sizes as low as 1, typical for streaming-dar settings, and is possible react and adapt changes in data statistics (concept drift/shift) without catastrophic forgetting. Major problems in such streaming settings are undesirable local optima during early training phases and numerical instabilities due to high data dimensionalities. We introduce an adaptive annealing procedure to address the first problem, whereas numerical instabilities to the standard GMM log- likelioodoood. Experiments on a variety of visual and non- visual benchmarks show that our SGD approach can be trained completely without, for instance, k-means based centroid initialization, and compares to a favorably (EM) stochastic EM (sEM), which it
eliminated by using exponential-free approximation
online variant of Sxspetation-Maxxiiaaaa outperforms by large margin for very high-dimensionall data.
INTRODUCTION
This contribution focuses Gaussian Mixture Models (GMMs), which represent probabilistic unsuper- vised model for clustering and density estimation and allowing sampling and outlier detection. GMMs have been used in wide range of scenarios, e.g., e.g., Melnykov & Maitra (2010). Commonly, free param- eters of GMM are estimated by using the Exppecaaioon-Maxinzzatti al., 1977).
(EMs) algorithm (Dempster) does not require learning rates and automatically enforces all GMM constraints. popular online variant is stochastic Expectation Maximization (sEM) (Cappé & Moulines) (2009). which can be trained mini-batch wise and is, thus, more suited for large datasets streaming data.
1.1 MOTIVATION
Intrinsically, EM is a batch-type algorithm. Memory requirements therefore become excessive for large datasets. In addition, streaming-date scenarios require data samples to be processed one by one, which is impossible for batch-type algorithm. Moreover, data statistics may be subject to changes over time (concept drift/shift), to which the GMM should adapt. In such scenarios, an online, mini-batch type of optimization such as SGD is attractive, it can process samples one by one, has modest, fixed memory requirements and can adapt to changing data statistics.
1.2 RELATED WORK
Online FM technique for performing EM mini-batch wise, allowing to process large datasets. One branch of previous research Newton et al. (1966): Lange (1995): Chen et al. (2018) has been devoted to the development of stochastic Expectation Maximization (sEM) algorithms that reduce the original EM method in the limit of large batch sizes. The variant of Cappé & Moulines (2009) widely used due to its simplicity and efficiency for large datasets. These approaches come at the price of additional syper-parametery (e.g., learning rate advantage of EM over SGD. Another common approach is to modify the EM algorithm itself by. e.g., including heuristics for adding, splitting and merging centroids Vlassis & Likas (002)) Engell
size). thus, removing key

---

Under review as a conference paper at ICLR 2021
& Heinen 2010); Pinto & Engel (2015): Cederborg al. (2010) Song & Wang (2005 Kristan et al. (2008): Vijayakumar et al. (2005). This allows GMM- like models to be trained by presenting one sample after another. The models work well in several application scenarios, but their learning dynamics are impossible to analyze mathematically. They also introduce a high number of parameters. Apart from these works, some authors avoid the issue of extensive datasets by determining smaller "core sets" of representative samples and performing vanilla EM Feldman et al. (201T).
SGD for training GMMs has, as far as we know, been recently treated only by Hosseini & Sra (2015) 2019). In this body of work, GMM constraint enforcement is ensured by using manifold optimization techniques and re- ce-paarmeeeeiizaiinnrggllaaaaaie are introduced. The issue of local optima is sidestepped by k-means type centroid initialization, and the used image datasets are low-dimensional (36 dimensions). Additionally, enforcing positive definiteness constraints by Cholesky decomposition is discussed.
Thereby, additional haperrparameters
Annealing and Approximation approaches for G.MMs were proposed by Verbeck et al. (2005);
Pinheiro & Bates (1995); Ormoneit & Tresp Dognin et al. (2009). However, the regularizers proposed by Verbeek et al. (2005); Ormmonei & Tresp (9988 significantly differ from our scheme.
GMM log-likelihood approximations, similar the used here, discussed in, e.g.. Pinheiro & Bates (1995) and Dognin et al. (2009), but only in combination with EM training.
GMM Training in ffigh-Dimensional Spaces is discussed in several publications: A conceptually very interesting procedure is proposed by Ge et al. (015)) It exploits the properties of high- dimensional spaces in order achieve learning with a number of samples that is polynomial in the number of Gaussian components. This is difficult to apply in streaming settings, since higher-order moments need to be estimated beforehand, and also because the number of samples usually cannot be controlled in practice. Training GMM-like lower-dimensional factor analysis models by SGD on high-dimensional image data is successfully demonstrated in Richardson & Weiss (2018), avoiding numerical issues, but, again, sidestepping the local optima issue by using k-means initialization. The numerical issues associated with loglliillhoodd computation in high- high-dimensional spaces are generally mitigated by using the "logsumexp" trick NNeessen & Sun (2016), which is, however, insufficient for ensuring numerical stability for particularly high-dimensional data, such as images.
1.3
GOALS AND CONTRIBUTIONS
The goals of this article are to establish GMM training by SGD as a simple and scalable alternative to sEM in streaming scenarios with potentially high- dimensiiiaii data. The main novel contributions are:
proposal for numerically stable GMM training by SGD that outperforms sEM for high data dimensionalities
an automatic annealing procedure that ensures SGD convergence from wide range of initial conditions without prior knowledge of the data (e.g., no k-means initialization) which is especially beneficial for streaming data computationally efficient method for enforcing all GMM constraints in SGD
Apart from these contents, we provide a publicly available TensorFlow implementation !
DATASETS
We use a variety of different image-based datasets as well as a non-image dataset for evaluation purposes. All datasets are normalized to the [0, 1] range.
MNIST (LeCun al. 1998) contains gray scale images, which depict handwritten digits from 0 to in resolution of 28 x28 pixels the common benchmark for computer vision systems.
SVHN (Wang et al.) 2012)) contains color images of house numbers (0-9, resolution 32 x 32). MashionMNiST Xiao et al. 2017) contains gray scale images of 10 cloth categories and is considered as more challenging classification task compared to MNIST.
Fruits 360 (Murean& & Oltean 2018) consists of colored pictures showing different types of fruits (100 x 100 3 pixels). The ten best-represented classes are selected from this dataset. Devanagari (Acharya et al. 2016) includes gray scale images of handwritten Devanagari letters with a resolution of 32 x32 pixels the first 10 classes are selected.
NotMNIST Yaroslav Bulatov 2011) is a gray scale image dataset (resolution 28 x 28 pixels) of
https//ggthhaccmmmmmmmmmlllll1lu
2

---

Under review as a conference paper at ICLR 2021
letters from A to extracted from different public available fonts.
ISOLET (Cole & Fanty] 1990) a non-image dataset containing 7 797 samples of spoken letters recorded from 150 subjects. Each sample encoded and is represented by 617 float values.
GAUSSIAN MIXTURE MODELS
GMMs probabilistic models that intend to explain the observed data X [z,} by expressing their density as a weighted mixture of Gaussian component densities N(;;; ph. Px) (=)) p(æ) = ΣE (x). We work with precision matrices P = Σ instead of covariances E. This realized by optimizing the (incomplete) log- oggllikelhhood
En
log Σ (xNN
3.1 GMMS AND SGD
GMMs require the mixture weights to be normalized: Σx T = and the precision matrices to be positive definite: x" P.æ 0 Va. These constraints must be explicitly enforced after each SGD step:
Weights Ty are adapted according to Hosseini Sra (2015), which replaces them by other free parameters & from which the TT are computed that normalization is ensured:
exp(Ex)
(2)
exp(Es)
Precision Matrices need to be positive-definite, so we re-parameterize these as P. p== PD= D&). where the upper-diagonal matrices D₁ result from Cholesky decomposition. Consequently, det Σg det P. (de Dx)) - (tr(Ds)) can be computed efficiently. To avoid re- computing the costly Cholesky decomposition of the Ps at every iteration, perform it on the initial precision matrices and just erase the elements below the diagonal in the D₃ after each gradient step.
3.2 MAX-COMPONENT APPROXIMATION FOR GMMS
The log- likelihoodd Eq. ( difficult to optimize by SGD (see Sec. 3.3 This is why we intend to find a lower bound that we can optimize instead. A simple scheme is given by
En log Σ (Z)))
> Ex log max; (TANN ( -Č= log (Te- (zw))
(3)
where A = arg max; (T.N)(( This is what we call the max-component approximation of Eq. (3). In contrast to the lower bound that is constructed for EM-type algorithms, this bound is usually tight. The advantages of Ć are the avoidance of local optima SGD, and the elimination of exponentials causing numerical instabilities for high data dimensions. The "logsumexp" trick is normally employed with GMMs rectify this by factoring out the largest component probability Nr. This mitigates, but does not avoid numerical problems when distances are high. To give an example: we normalize the component probability M₁. = (using 32-bit floats) by the highest probability Ng. 3 and obtain N 104 which is numerically problematic.
3.3 UNDESIRABLE LOCAL OPTIMA IN SGD TRAINING
An issue when performing SGD without k-means initialization concerns undesirable local optima. Degenerate Solutions occur when naively optimizing L by SGD (see Fig. Tap. All components have case all gradients vanish (see App. A.3 proof). These solutions subset of components is updated by SGD, thereby breaking the symmetry between components.
the same weight, centroid and covariance matrix: मे R k. #4 = E Σ (xv(X) ) A k, in which
avoided by Ł, since only
Single/Sparse-Compont Solutions occur when optimizing Ć by SGD (see Fig. lb). They are characterized by or several components (k:) that have large weights with centroid and precision
3

---

Under review as a conference paper at ICLR 2021
matrices given by the mean and covariance of significant subset X X of the data X: TT, > μk, = X Es; = Cov (X,, whereas the remaining components & are characterized by T Re 0. Mg = -== Pa Pt 0). Thus, these unconverged components are almost never best-matching components k" The max operrttiin in causes gradients like ac to contain ŌÆ (see App. A.3 tffe
This implies that they are non-zero only for the best-matching component k*. Thus the gradients of unconverged components vanish, implying that they remain their unconverged state.
0.94A
0.40
0,X
01000
0,010
0.110
0.010
0.010
0.00
0.034
9.001
0.500
0.000
0,819
0,024
0,00L
0.024
0.00
3.001
0.040
0,640
0.640
0,0%4
0,0%
DUNG
0.440
0.000
0.000
0.00
00X
a) degenerate solution (optimizing L)
b) component solution (optimizing 2)
Figure 1: Undesirable solutions during SGD, visualized for MNIST with component weights Tk-
3.4 ANNEALING PROCEDURE FOR AVOIDING LOCAL OPTIMA
Our approach for avoiding undesirable solutions is to punish their characteristic response patterns by a modification of É, the smoothed max-component log-likelihood L
L" E,,
Σ (o) (1olgg(TT,,, (T.NN((E.)) = En Σ gh (() (T)N)
(4)
The entries of the are computed by a Gaussian function centered on component h with common spatial standard deviation where we assume that the K components are arranged on v /K' x VK' grid with 2D Euclidean metric (see App. [A.4p. Eq. essentially represents a smoothing of the log (x=) with 2D convolution filter (we use periodic boundary conditions). Thus, Eq. (4) maximized if the log probabilities follow uni-modal Gaussian profile of spatial variance o which heavily punishes single-componee solutions that have locally delta-like response.
Annealing starts with large values of (l) = and reduces over time to an asymptotic small value of 0 000, thus, smoothly transitioning from _ in Eq. (4 into P Eq. B).
Annealing Control is ensured by adjusting a, which defines an effective upper bound on 2 (see App. A.2 for proof). This implies that the loss will be stationary once this bound reached, which we consider a suitable indicator for reducing a. We implement an annealing control that sets 0.97 whenever the loss considered sufficiently stationary. Stationarity is detected by maintaining an exponentially smoothed average E(t) = (1 a)e(t - 1) + aL on time scale a.
Every 1 iterations, we compute the fractional increase of L" as
((t) - f(t - a -1 E(t - -1 (t and consider the loss stationary iff A < 8 (the latter being free parameter). The choice of the time constant for smoothing La is outlined in the following section.
3.5 TRAINING PROCEDURE FOR SGD
Training GMMs with SGD is performed by maximizing the smoothed max-component log-likelihood C" from Eq. (4). At the same time, enforce the constraints the component weights and covariances described in Sec. 3.1 and transition from 2 into L by annealing (see Sec. 3.4 SGD requires a learning rate e to be set, which in turn determines the parameter a (see Sec. 3.4) as a = since stationarity detection should operate time scale similar to that of SGD. Cholesky matrices
4

---

Under review as a conference paper at ICLR 2021
D3 are initialized to max and are clipped after each iteration so that diagonal entries are in the range [0, D This is necessary to avoid excessive growth of precisions for data entries with vanishing variance, e.g., pixels that are always black. Weights are uniformly initialized to = K centroids the range + (see Alg. 1 for summary). Please note that our SGD approach requires centroid initialization by k- k-means, it is usually recommended when training GMMs with EMM We discuss and summarize good practices for choosing aper-parameters Sec.
Algorithm 1: Steps of SGD-GMM training.
Data: initializer values: n², K, eoltoo, aolo ô and data X Result: trained GMM model
μ f + 20-dd(''' --"=-="''')
f- 1/K, P P+-IDmus mDxx
ou. ct 6o
initialization
2 forall < do
training loop A.4
g(t) f preate-annealing mask(o,t)
see Sec.
3.4 and App.
pu(t), P(t), (t) ač
E aco By + pu(-1),, e "P + P(t-1),
0t + x(t-1)
I/ SGD update Sec. // see Eq.
P(t) precisions clipppinggnippppp(PP... Dmax)
3.5
<(t) + normalization(((())
L (1)+ (1 (1-a))al) a)](t-1) + (x(())
I/ sliding log-likelihood Sec. see Eq.
annealing update iteration then
A < then (t) <- 0.90(t - 1), ((t) + 0.9
3.6 TRAINING PROCEDURE FOR STOCHASTIC EXPECTATION MAXIMIZATION
We use sEM as proposed by Cappé & Moulines (0099. We choose the step size of the form P: = (t 1)-0.5+a with E [0, 0.5), Po and enforce p(t) 2 Poo- Values for these parame- ters are determined via a grid search in the ranges po E (0.01, 0.05, 0.1}, a (0.01, 0.25, 0.5) and Poo {0.01, 0.001, 0.0001 Each sEM iteration uses batch size B. Initial accumulation of suf- ficient statics is conducted for 10% of an epoch, but not when re-training with new data statistics.
Parameter initialization and clipping of precisions is performed just for SGD, see Sec. 3.5
3.7 COMPARING SGD AND SEM
Since sEM optimizes the log-likeli L, whereas SGD optimizes the annealed approximation L the comparison of these measures should be considered carefully. We claim that the comparison is fair assuming that i) SGD annealing has converged and ii) GMM responsibilities are sharply peaked so that single component has responsibility of x 1. It follows from i) that L x É and ii) implies that e x L. Condition ii) is usually satisfied to high precision especially for high-dimensional inputs: if it is not, the comparison is biased in favor of sEM, since >> L by definition.
EXPERIMENTS
Unless stated otherwise, the experiments in this section will be conducted with the following parameter values for sEM and SGD (where applicable): mini-batch size B = 1, K== x 8, μ² =²= 0.1,, op = x 0.01, 0.001, D 20. Each experiment is repeated 10 times with identical parameters but different random seeds for parameter initialization. See Sec. | justification of these choices. Due to input dimensionality, all precision matrices must be taken to be diagonal. Training/test data are taken from the datasets mentioned in Sec.
4.1 ROBUSTNESS OF SGD TO INITIAL CONDITIONS
Here, we train GMMs for three epochs on classes 9 for each dataset. We use different random and non-random initializations of the centroids and compare the final log-likelihood values. Random centroid initializations are parameterized by μr E {0.1, 0.3, 0.5}, whereas non-random initializations class 0 (one epoch). The latter is done to
defined by centroids from previous training have a non-random centroid initialization that as dissimilar as possible from the training data. The initialization of the precisions cannot be varied because empirical data shows that training converges
5

---

Under review as a conference paper at ICLR 2021
to undesirable solutions if the precisions not initialized to large values. While this will have to be investigated further, we find that convergence to near eear-idettiaal levels, regardless of centroid initialization for all datasets (see Tab. for more details).
4.2 ADDED VALUE OF ANNEALING
To demonstrate the beneficial effects of annealing, we perform experiments on all datasets with annealing turned off. This is achieved by setting op = dxo. This invariably produces sparse- component solutions with strongly inferior bog-liielihoods after training, please refer to Tab. [
Table I: Effect of different random and non-random centroid initializations of SGD training. Given are the means and standard deviations of final log- oggliiklssoods (10 repetitions per experiment). To show the added value of annealing, the right-most column indicates the final log-likelihoods when annealing is turned off. This value should be co.pared to the leftmost entry in each row where annealing is turned on. Standard deviations in this case where very small so they are omitted.
Initialization
random
non-random
no ann.
/ = 0.1
J² =
0.3
14
0.5
init class 0
/ = 0.1
Dataset
mean
std
mean
std
mean
atd
mean
std
meon
MNIST
205.47
1.08
205.46
0.77
205.68
0.78
205.37
0.68
124.1
FaSHionMNIST
231.22
1.53
231.58
2.84
231.00
1.11
229.59
0.59
183.0
NotMNIST
- -18.41
1.77
- 48.59
1.56
-48.32
1.13
-49.37
2.32
-203.8
Devanagari
15.95
1.59
15.76
1.34
-17.01
1.11
-22.07
4.59
-263.4
Fruits 360
095.80
98.02
00000.70
127.00
036.25
122.06
10912.79
1 1777661
331.2
SVHN
1328.06
0.94
1327.99
1.59
1 1328.40
1.17
1 1327.80 80
0.94
863.2
ISOLET
354.34
0.04
354.36
0.04
354.36
0.04
354.20
0.05
201.5
4.3 CLUSTERING PERFORMANCE EVALUATION
To compare the clustering performance of sEM and GMM the Davies-Bouldis score (Davies & Bouldin 1979) and the Dunn index (Dunn 1973) are determined. We evaluate the grid-search results to find the best parameter setup for each metric for comparison. sEM is initialized by k-means to show that our approach does not depend on parameter initialization. Tab. 2 indicaties that SGD can egalize sEM performance (see also App. A.5).
Table 2:
Clustering performance comparison of SGD and sEM training using Davies Bouuldn score (less is better) and Dunn index (more is better). Results are in bold face whenever they are better by more than half a standard deviation.
Metric Algo.
Davies-Bouldin score sEM
Dunn index
SGD
SGD
sEM
Dataset
meon
std
meon
std
std
mecn
std
MNIST
2.50
0.04
2.47
0.04
0.18
0.02
0.16
0.02
FashionMNIST
2.06
0.05
2.20
0.04
0.20
0.03
0.19
0.02
NotMNIST
2.30
0.03
2.12
0.03
0.15
0.03
0.14
0.04
Devanagari
2.60
0.04
2.64
0.02
0.33
0.01
0.27
0.04
SVHN
2.34
0.04
2.41
0.03
0.15
0.02
0.15
0.02
4.4 STREAMING-DAA DATA EXPERIMENTS WITH CONSTANT DATA STATISTICS
We train GMMs for three epochs (enough for convergence all cases) using SGD and sEM on all datasets as described in Secs. 3.5 and 3.6 The resulting centroids of our SGD-based approach are shown in Fig. [ whereas the final loss values for SGD and sEM are compared in Tab. 3 The centroids for both approaches are visually similar, except for the topological organization due to annealing for SGD, and the fact that in most experiments, some components do not converge for sEM while the others do. Tab. 3 indiccaes that SGD achieves performances superior sEM in the majority of cases, in particular for the highest-dimensiona. dimensional datasets (SVHN: 3 072 dimensions and Fruits 360: 30 000 dimensions).

---

Under review as a conference paper at ICLR 2021
Table Comparison of SGD and sEM training all datasets streaming- scenario. Each time mean log-likeline liglliies (10 repetitions) at the end of training, and their standard deviations are presented. Results are in bold face whenever they are higher by more than half a standard deviation. Additionally, the averaged maximum responsibilities (pk-) for test data are given for justifying the max-component approximation.
Algorithm
SGD
SEM
Dataset
Ø max Pi:-
mean
std
mean
std
MNIST
0.992 674
216.60
0.31
216.82
1.38
FaShionMNIST
0.997 609
234.59
2.28
222.95
6.03
NoMMNIST
0.998 713
34.76
1.16
- 40.09
8.90
Devanagari
0.999 253
14.65
1.09
13.46
6.16
Fruits 360
0.999 746
11754.32
75.63
5 4-43.00
201.60
SVHN
0.998 148
1329.80 83
0.80
176.007
16.91
ISOLET
0.994 069
354.21
0.01
354.55
0.37
MNIST
Devanagari NotMNIST Fruits 360 Figure 2: Exemplary results for centroids learned by SGD, trained full images.
SVHN
MNSST
Visualization of High-dimensional sEM Outcomes Fig. 3 was obtained after training GMMs by sEM on both the Fruits 360 and the SVHN dataset. It should be compared to Fig. L where an identical procedure was used to visualize centroids of SGD-trained GMMs. It is notable that the effect of unconverged components does not occur at all for our SGD approach, which is due to the annealing mechanism that "drags" unconverged components along.
Figure 3: Visualization of centroids after exemplary training runs (3 epochs) on high -igh--imnnsional datasets for sEM: Fruits 360 (left, 30 000 dimensions) and SVHN (right, 000 dimensions). Com- ponent entries displayed "as is", meaning that low brightness means low RGB values. Visibly, many GMM components remain unconverged, which is analogous to a sparse-component solution and explains the low log-likelihood values especially for these high-dimensional datasets.
DISCUSSION AND CONCLUSION
The Relevance of this Article is outlined by the fact that training GMMs by SGD was recently investigated in the community by |Hosseini & Sra (2015) 2019). We go beyond, since our approach does not rely on off-line data-driven model initialization, and works for high-dimensional streaming dhea.The presented SGD scheme is simple and very robust to initial conditions due to the proposed annealing procedure, see Sec. 4.1 and Sec. 4.2 In addition, our SGD approach compares favorably to the reference model for online EM (Cappé & Moulines 2009) in terms of achieved log-likelihode likelihoods, which was verified on multiple real-world datasets. Superior SGD performance observed for the high-dimensional datasets.
7

---

Under review as a conference paper at ICLR 2021
Analysis of Results suggests that SGD performs better than sEM average, see Sec. 4.4 although the differences are very modest. It should be stated clearly that it cannot be expected, and is not the goal of this article, to outperfom sEM by SGD in the general case, only to achieve similar performance. However, if sEM is used without, e.g., k-means initialization, components may not converge (see Fig. for a visual impression) for very high-dimensional data like Fruits 360 and SVHN datasets, which is why SGD outperforms sEM in this case. Another important advantage of SGD over sEM is the fact that no grid search for finding certain whereas sEM has a complex and unintuitive dependency on pu, Poo and Cug
values is necessary,
Small Batch Sizes and Streaming Data are possible with the SGD-based approach. Throughout the experiments, we used a batch size of 1, which allows streaming- data processing without the need to store any samples at all. Larger batch sizes are, of course, possible and increase execution speed. In the experiments conducted here, SGD (and sEM) usually converged within the first two epochs, which is a substantial advantage whenever huge sets of data have to be processed.
No Assumptions About Data Generation are made by SGD in contrast to the EM and sEM algorithms. The latter guarantee that the loss will not decrease due to an M-step. This, however, assumes non-trivial dependency of the data proof). In contrast, SGD makes no such tard-to-verify assumptions, which is a rather philosophical point, but may be an advantage in certain situations where data are strongly non-Gaussian.
unobservable latent variable (see A for
Numerical Stability is assured by our SGD training approach. does not optimize the log-likelhhood but its max-component approximation. This approximation contains no exponentials at all and is very well justified by the results of Tab. 3 which show that component probabilities are very strongly peaked. In fact, is the gradient computations where numerical problems (e.g.. NaN values) occurred. The "logsumexp" trick mitigates the problem, but does not eliminate it (see Sec. 3.2). It cannot be used when gradients are computed automatically (what most machine learning frameworks do).
Hyper-Parameter Selection Guidelines are as follows: the learning rate e must be set by cross- validation (a good value is 0.001). We empirically found that initializing precisions to the cut-off value Dmaax and an uniform initialization of the are beneficial, and that centroids are best initialized to small random values. A value of Dmmax = 20 always worked in the experiments. Generally, the cut-off must be much larger than the inverse of the data variance. In many cases, it should be possible to estimate this roughly, even in streaming settings, especially when samples are normalized. For density estimation, choosing higher values for K leads to higher final log- ligllikelihoods (validated in App. App.AÃ.6) [A. For clustering, K should be selected using standard techniques for GMMs. The parameter ō controls loss stationarity detection for the annealing procedure and was shown to perform well for 5 = 0.05. Larger values will lead to faster decrease of (4), which may impair convergence. Smaller values are always admissible but lead to longer convergence times. The annealing time constant a should be set to the GMM learning rate or lower. Smaller values of a lead to longer convergence times since (t) will be updated less often. The initial value dg needs to be large in order to enforce convergence for all components. A typical value 0. 0.25v/ K The lower bound on a. doc should be as small as possible to achieve high log-likelihoods (e.g.. 0.01, see proof in App. A.2
6 OUTLOOK
The presented work can be extended in several ways: First of all, annealing control could be simplified further by inferring good 5 values from a. Likewise, increases of a might be performed automatically when the loss rises sharply, indicating a task boundary. As we found that GMM convergence times grow linear with the number of components, we will investigate hierarchical GMM models that operate like Convolutional Neural Network (CNN), and in which individual GMMs only see local patch of the input and can therefore have low K. Lastly, we will investigate a replay by SGD- trained GMMs for continual learning architectures. GMMs could compare favorably to Generative Adversarial Nets (Goodfellow et al. 2014) due to faster training and the fact that sample generation capacity can be monitored via the logllikeli
8

---

Under review as a conference paper at ICLR 2021
REFERENCES
Shailesh Acharya, Ashok Kumar Pant, and Prashnna Kumar Gyawali. Deep learning based large scale handwritten Devanagari character recognition. SKIMA 2015- 9th International Conference on Software, Knowledge, Information Management and Applications, 2016. doi: 10.1109/SKIMA 2015.7400041.
Olivier Cappé and Eric Moulines. On-line expectation expectatio-maxxiiizttooi algorithm for latent data models. Journal of the Royal Statistical Society, Series B: Statistical Methodology, 71(3):593-6133 2009.
ISSN 13697412. doi: 10. (0.1111/: 10111//.11467 1467-9868 2009,006988
Thomas Cederborg, Ming Li, Adrien Baranes, and Pierre Yves Oudeyer. Incremental local online Gaussian Mixture Regression for imitation learning of multiple tasks. IEEE/RSJ 2010 International Conference on Intelligent Robots and Systems, IROS 2010 Conference Proceedings, pp. 267-274, 2010. doi: 10. 10.1100 0009//ROS.2000.6655
Jianfei Chen, Jun Zhu, Yee Whye Teh, and Tong Zhang. Stochastic expectation maximization with variance reduction. Advances in Neural Information Processing Systems, 8-Deembbr (NeurIPS),7967-797,, 2018.
2018. ISSN 10495258.
Ronald Cole and Mark Fanty, Spoken letter recognition, pp. 385-390, 1990. doi: 10.31 1..3115/116580 116725.
David L. Davies and Donald W. Bouldin. A Cluster Separation Measure. IEEE Transactions on Pattern Analysis and Machine Intelligence, PAMI-1(2:224-227, 1979. ISSN 01628828. doi: 10.1109/TAMMII1993 19T1AAA99999977669999
A. P. Dempster, N. M. Laird, and D. B. Rubin. Maximum Likelihood from Incomplete Data Via the EM Algorithm. Journal of the Royal Statistical Society: Series B (Methodological), 39(1):1-22, 1977. doi: 10.1111/j22517-61
1977.100111 ib016600.
Pierre L. Dognin, Vaibhava Goel, John R. Hershey, and Peder A. Olsen. A fast, accurate approximation to log likelihood of Gaussian mixture models. ICASSP, IEEE International Conference on Acoustics, Speech and Signal Processing- Proceedings, (3):3817-3820, 2009. ISSN 15206149. doi: 10.1 1091 NASSSP:2009,4960453
J. C. Dunn. A fuzzy relative of the ISODATA process and its in detecting compact well- separated clusters. Journal of Cybernetics, 3(3):32-57, 1973. ISSN 00220280. doi: 1080/ 0196972730854604.
Paulo Martins Engel and Milton Roberto Heinen. Incremental learning of multivariate Gaussian mixture models. Lecture Notes in Computer Science (including subseries Lecture Notes in Artificial Intelligence and Lecture Notes Bioinformatics), 6404 LNA1:82 91, 2010. ISSN 03029743. doi:
10.1007/978-3-6422 1613--4 9.
Dan Feldman, Matthew Faulkner, and Andreas Krause. Scalable training of mixture models via coresets. Advances in Neural Information Processing Systems 24: 25th Annual Conference on Neural Information Processing Systems 2011, NIPS 2011, pp. 1-9, 2011.
Rong Ge, Qingqing Huang, and Sham M. Kakade. Learning mixtures of gaussians in high dimensions. Proceedings of the Annual ACM Symposiuum on Theory of Computing, 14-17-lum-200157761-770. 2015. ISSN 07378017. doi: (0.1145/2744639.2746611
Ian Goodfellow, Jean Pougct-Abadic, Mehdi Mirza, Bing Xu, David Warde-Farley, Sherjil Ozair, Aaron Courville, and Yoshua Bengio. Generative adversarial nets. Z. Ghahramani, M. Welling, C. Cortes, N. D. Lawrence, and K. Q Weinberger (eds.), Advances in Neural Information Process- ing Systems 27, pp. 2672-26800 Curran Associates, Inc., 2014, URL httt:://papee: nips ce/ cc/ htlp:/aprr//// 542 httprr/eeeeeeeeerreeeeeeeeeeeeee
pdf
Reshad Hosseini and Suvrit Sra. Matrix manifold optimization for Gaussian mixtures. Advances in Neural Information Processing Systems, 2015 111.5-JJaaar::100-11 10-9188, 2015. ISSN 10495258.
Reshad Hosseini and Suvrit Sra. An alternative to em for gaussian mixture models: Batch and stochastic riemannian optimization. Mathematical Programming, pp. 2019.

---

Under review as a conference paper at ICLR 2021
Matej Kristan, Danijel Skocaj, and Aleš Leonardis. Incremental learning with gaussian mixture models. In Computer Vision Winter Workshop, pp. 25-32, 2008
Kenneth Lange. A Gradient Algorithm Locally Equivalent to the EM Algorithm. Journal of the Royal Statistical Society: Series B (Methodological), 57(2)422-437, 1995. ISSN 0035-9246. doi:
Yann LeCun, Léon Bottou, Yoshua Bengio, and Patrick Haffner. Gradient-based learning applied to document recognition. Proceedings of the IEEE, 86(11):2278 2323, 1998. ISSN 00189219. doi:
10.1109/5.72671. /tt/tttp://eeeeelorreeeee wwwpppooreeee.ee..eo.e.eeeeeeeee htt.://e..rroococmmet////22353
Volodymyr Melnykov and Ranjan Maitra. Finite mixture models and model-based clustering. Statistics Surveys, 4(00etobe):80-116 2010. ISSN 19357516. doi: 10.1214/09-SS05,
Horea Murean and Mihai Oltean. Fruit recognition from images using deep learning. Acta Universi-
tatis Sapientiae, Informatica, 10(1):26-42) 2018. ISSN 2066-7760, doi: 10.2478/ausi-2018- 2008-002..
URL URLhttp: /tap:v//axxi...m org/ abs 0elpo//112.005880 00580
J. Newton, D. M. Titterington, A. F. M. Smith, and U. E. Makov. Statistical Analysis of Finite Mixture Distributions. Biometrics, 42(3):679, 1986. ISSN 0006341X. doi: 10.2307/2531224.
Frank Nielsen and Ke Sun. Guaranteed bounds on information coommation-theortticc measures of univariate mixtures using piecewise log-sum-ext inequalities. Entropy, 18(12), 2016. ISSN 10994300. doi: 10.3390/e 81 20442.
Dirk Ormoneit and Volker Tresp. Averaging, maximum penalized likelihood and Bayesian estimation for improving Gaussian mixture probability density estimates. IEEE Transactions on Neural Networks, (())6339-65) 1998. ISSN 10459227. doi: 10.1 10.1109/72.701177.
Jose C. Pinheiro and Douglas M. Bates. Approximations to the log-likelihood function in the nonlinear mixed-effects model. Journal of Computational and Graphical Statistics, 1(1):1 1995. ISSN 15372715. doi: 10.1080/106 002111616600 19900.0199.11044666
Rafael Coimbra Pinto and Paulo Martins Engel. A fast incremental Gaussian mixture model. PLoS ONE, 10(10), 2015. ISSN 19326203. doi: 10. 13.137 //11////ouaall.poee.... 0139911....91......
Eitan Richardson and Yair Weiss. On GANs and GMMs. Advances in Neural Information Processing Systems, 2018- December(NeurIPSS5554
2018. ISSN 10495258.
Mingzhou Song and Hongbin Wang. Highly efficient incremental estimation of Gaussian mixture models for online data stream clustering. Intelligent Computing: Theory and Applications III. doi: 10. 10.1117/12.601722.
5803:174, 2005.
J. Verbeek, N. Vlassis, and B. J.A. Kröse. Self-organizing mixture models. Neurocomputing, 63 (SPEC. 18S.):99-123, 2005. ISSN 09252312. doi: 10.10166j.neucom.2000 2004.00..0000000000..
Sethu Vijayakumar, Aaron D' D'Souza, and Stefan Schaal. Incremental online learning in high di- mensions. Neural Computation, (712)):602-26344 2005. ISSN 08997667. doi: 10.1162/
089976605774320557.
Nikos Vlassis and Aristidis Likas. A greedy EM algorithm for Gaussian mixture learning. Neural Processing Letters, 15(1):77-87, 2002. ISSN 13704621. doi: 10.1023/4:011338811113
Tao Wang, David J Wu, Adam Coates, and Andrew Y Ng. End-to-end text recognition with convolutional neural networks. In Proceedings of the 21st International Conference on Pattern Recognition (ICPR2012), pp. 3304-3308. IEEE, 2012.
Han Xiao, Kashif Rasul, and Roland Vollgraf. Rashion-MNISTI a Novel Image Dataset for Bench- marking Machine Learning Algorithms. pp. 1-6, 2017. URL htttp://arxiw.orr////3
1708.07741
Yaroslav Bulatov. Machine Learning, etc notMNIST dataset, 2011. URL htLntt:pp/yaarlaava bl www.spo.t....mm....o....o.o..... com br. http://t...oooooo///////3/ attp html
10

---

Under review as a conference paper at ICLR 2021
A SUPPLEMENTARY MATERIAL
A.1 ASSUMPTIONS MADE BY EM AND SGD
The EM algorithm assumes that the observed data samples {œ,, depend on unobserved latent variables 2n in a non-trivial fashion. This assumption is formalized for a GMM with K components by formulating the complete-data likelihood in which 2n 2n {O...., K 1) is scalar:
(6) where we have defined Nig (ax,) = (x 0% p2) for brevity. is assumed that the 2n are unob- servable random variables whose distribution is unknown. Marginalizing them out gives us the
p( 2n) N2. (2:
incomplete-data likelihood p(æ,,) Σx p(In, Zn). The derivation of the EM algorithm starts out with the total incomplete-dan data log- loellikel
log p(X log II [p(xn)
Σ log p(x)))
Σ log Σ plans on = k))
log Σ p(an 1:) p( k)
p(zn -=n=l
Due to the assumption that C is obtained by marginalizing out the latent variables, an explicit dependency on 2n can be re-introduced. For the last expression, Jensen' inequality can be used to construct lower bound:
Σ log Σ plzz = k)) P( Zn = k) 12
L =~
p(zn k)
k
(8)
p(z, = log p( Zn = k)
plzn k)
Since the realizations of the latent variables are unknown, we can assume any form for their distribu- tion. In particular, for the choice p(2n) ~ p(æ,, 2n), the lower bound becomes tight. Simple algebra and the fact that the distribution p(za) must be normalized gives us:
p(zn- k, k,)
p(zn = k))
p(xn.)
p(zn k3 (2)))))))
plzn n ki 2kn)
(9)
P(z = 1, x.n)
(2
Σ, (NN)
where we have used Eq. (6) in the last step. p(2n A k/lr,) is quantity that be computed from data with no reference to the latent variables. For GMMs it is usually termed responsibility and we write it p(=,. k| )..) 7nk.
However, the construction of a tight lower bound, which is actually different from L, only works when p(æn, 2n) depends non-trivially on the latent variable 2n. If this is not the case, have p(Cn, 2n) = K -))((x,) and the derivation of Eq. (8) goes down very differently:
p(un: 2n k)
C
log logp(w))) 2 Σ Σ
ΣΣ
p(zn = log
p(zn = k)
K
K- (<1p(()))))))))
Σ p(zn = k)) log
plzn
(10)
Σ log K (((ppa))))
Σ p(zn = k) logp(zn k)
Σ (logp(x.)) (log K H((zl)
11

---

Under review as a conference paper at ICLR 2021
where H. represents the Shannon entropy of p(z). The highest value this can have is log K for an uniform distribution the 2n, finally leading to a lower bound for L of
C 2 Σ (log((z))
(11)
which is indeed tight, but trivial, and thus does not simplify the problem at all. In particular, no closed-form solutions the associated extreme value problem can computed for this case.
This shows that optimizing GMMs by Expectation-Maxiin -aximizzatinn assumes that each sample has been drawn from single element in set of K uni-modal Gaussian distributions. Which distribution is selected for sampling depends on latent random variable. On the other hand, optimization by SGD uses the incomplete-data log-likelihood as basis for optimization, without probabilistic interpretation. This may be advantageous for problems where the assumption of Gaussianity is badly violated, although empirical studies indicate that optimization by EM works very well a wide range of scenarios.
A.2 PROOF THAT DEFINES AN UPPER BOUND ON L"
Let us assume that SGD optimization has reached stationary point where the derivative w.r.t. all GMM parameters O. In this situation, we claim that the only way to increase the loss is by ać" manipulating a. We show here that do < <va > O, and that _ć de =0 for =0. This means that the loss can be increased by decreasing a, up the point where 0
For each sample, the 2D profile of log(m = ʃ is assumed to be radially symmetric, cen- tered on the best-matching component K- and decreases with the distance as a function of ||k k-ll- We thus have fa = fx(r) with r= || - k*1l. Passing to the continuous domain, the indices in the Gaussian "smoothing filter" Λ become continuous variables, and we have I + 9(lk - k" 11, o) g(r, or. Similarly, fa(r) + f(r). Using 2D polar coordinates, the smoothed max-component likelihood L becomes polar integral around the position of the best- matching component: L" ~ Juz (
We interested in the change of L when undergoes an infinitesimal change. It is trivial to show that for the special case of a constant log- probability profile, i.e., (r) = L, L does not depend a because Gaussians are normalized, and that the derivative w.r.t. a vanishes:
dL
l dr 12
72
02 -1) exp((
202 )
do
0
L I dr
12
,2
122 (- 02 -1) exp(-
2²
(12)
0² -1)
exp(
L
202
20
J0
LN - LP
where we have split the integral into parts where the derivative w.r.t. a is negative(N) and positive
(P). We know that / = P since the derivative must be zero for constant function f(r) L. due to the fact that Gaussians are normalized to the same value regardless of a.
For function fr) that satisfies f(r) > Ltr [ and f(r) Lvr Ela, oo, the inner and outer parts of the integral behave as follows:
2 g(r) -1) f g(r) (( 2 0
Ñ
( -1) LN
(13)
< -1) g(r) 2
2 g(r) ( -1) f() <
P
--1- LP
since f(r) is minorized//aaoorzedd by by assumption, and the contributions in both integrals have the same sign for the whole domain of integration. Thus, it is shown that, for o > 0
dc. - Ñ - LN LP =
(14)
do
and, furthermore, that this derivative is zero for σ = because La no longer depends on o for this case.
12

---

Thus, a indeed defines an upper bound L which can be increased by decreasing o. The assumption of pog-prbbabilities that decrease around the best matching unit (BMU) is reasonable since such profile maximizes É" All functions fr) that, e.g., decrease monotonically around the BMU, fulfill this criterion, where the precise form the decrease is irrelevant. This proof works identically when Under review as a conference paper at ICLR 2021
Taking everything into consideration, in situation where the log-likelihood La has reached stationary point for given value of a, we have shown that:
The value f" depends on a.
Without changing the log-probabilities, we can increase Ωe by reducing assuming that the log-probabilities mildly decreasing around the BMU.
Increasing L works long
At 0 the derivative becomes 0.
not resorting to integrals but using discrete sums.
A.3 LOG-LIKELIIOOD GRADIENTS
The gradients of read
σC
E, [P& (I, (4x) 6RA-]
0p
En [((P(( (0n 14:) -:1))) 6k
))
(15)
OP
OC
Tkk En, [ŌRA-].
0Tk
The gradients of L are obtained by replacing 6&g- by the standard GMM responsibilities "int For the case of degenerate solution when optimizing L. only single component k* has a weight close to 1, with its centroid and covariance matrix are given by the mean and covariance of the data: T& Mg. E[ P-- Cov(X). In this case, the gradients w.r.t μe and P vanish. The gradient w.r.t. Thk does not vanish, but is 6gs which vanishes after enforcing the normalization constraint.
13

---

Under review as a conference paper at ICLR 2021
A.4 VISUALIZATION OF 2D ANNEALING GRID
In Figure 4, three different states of the Oh are visualized, depending on o(l). Darker pixels indicate larger values. Each Ji is assigned to single GMM component k; which why the gh arranged on the same vK" x v K" grid place the components themselves. An intuitive interpretation of particular 9s. is that it encodes the contribution of neighbouring component log log-probabilities on the tog-probability of component & entering into nax-computation. Over time, (t) reduced (middle and right pictures) and thus only component & contributes. Please note that the grid we place the components on is periodic for simplicity, so the g& are themselves periodic.
Figure 4: Visualization of different centered Gaussian function gt states controlled by
14

---

Under review as a conference paper at ICLR 2021
A.5
DIAGRAMS OF CLUSTERING COMPARISON
MNIST
davies boalin soure
GMM
0.20
ARM
duns indlx
3.0
0.15
CIMM
2.5
0.10
KFM
50.000
100 00D
150.000
50 000
100 000
150 000
Iteration
Seration
FashionMNIST
daavi bonidnn anxe
3.0
+ GMM
duan inlox
0.15
2.5
0.10
GMM
SEM
0.05
000
0000 000
1 500 000
500 000
1 BOG EX10
500 000
iteration
steration:
NotMNIST
davins boalli soure
3.40
GMM
WFM
dann index
0.15
2.5
CIMM
0.05
000000
00000000
5000
000
500
1000
10000000
iteaatin
Beration
Devanagari
davina boaldi anoe
+ GIMM
CIMM
daana lnslx
0..
0.2
10000
200000
30.000
40.000
50.000
10000
20000
30000
40000 500000
itrration
Reration
SVHN
davieo bonbin seate
+ GMD
duna indox
3
+ GMM
30
100000
150
(1X
500000
100000
150
200
Iteration
steratoon
Figure 5: Trend of clustering capabilities for sEM and SGD trained GMM. Comparison by Davies- Bouldin score (smaller better) and Dunn index (higher better) for the datasets MNIST, Fash- ionMNIST, NotMNIST, Devanagari and SVHN. The lines visualize the average metric score/index values of 10 repetitions with its standard deviation.
15

---

Under review as a conference paper at ICLR 2021
A.6 EFFECT OF THE GMM COMPONENT NUMBER K
The number of Gaussian components K a key parameter for any GMM and has huge impact on performance. It should be stated clearly that this discussion is not specific any particular fashion of performing GMM training, be it by SGD, sEM or EM since all of these methods optimize the log-likelihood. This is why we do not propose particular way of choosing K for SGD-trained GMMs. For clustering, K can be chosen using standard techniques like BIC or AIC, in addition to priors depending on the data and the concrete application generally assumed that K should be set as high as possible since the real distribution of the data can be approximated in better detail. This set of experiments aims at showing empirically, for completeness, that this "bigger is better" relation holds when training GMMs by SGD. We use the hyper-parameter settings stated in Sec. 4 vary K and record the final log-likelihoods. Results are shown in Tab, and suggest very clear relationship between K and the (test) log-likelinee lieelikiihoodd obtained the end of training. We run the experiments for larger number of epochs to exclude that (non- >covvergence effects could play role here.
mind. For density estimation, is
Table 4: Log-likelihood values after 50 training epochs using different number of Gaussian compo- nents K E {25, 36, 49, 64, 813. Best values of the final fog-likelihood are in bold face.
Components
SGD K = 7
K=
K =6
K =
K 9
Dataset
MNIST FNSSi
205.4
215.2
222.7
229.6
215.9
250.5
266.4
273.9
281.8 8
Devanagari
15.7
0.7
20.9
21.3
66.7 7
Fruits 360
12.95.7
13003.0
12352.2
12 866.6
13193.7
16


Reviews:
[{"review_id": "m6aeLE-0Be8", "reviewer": "AnonReviewer1", "summary": "", "questions": "", "limitations": "", "rating": 5, "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": "It is clear that efficiectly training Gaussian Mixture Models with SGD on non-stationary, high-dimensional streaming data is very important for practical applications. The aim of this paper is quite good.  In fact, this paper proposes an annealing mechanism for the SGD algorithm and makes the experiments to compare the proposed training algorithm with sEM algorithm on several real-world datasets.  However, I have following major concerns: \n(1). The proposed training scheme does not require data-driven parameter initialization (e.g., k-means) , but the data-driven parameter initialization can improve the efficiency. So, I cannot consider this is an advantage of the proposed scheme. \n(2).  The  proposed annealing scheme is straightforward，and there is no deep analysis on its performance.\n(3).  According to the experimental results, I cannot find out that the proposed training scheme is remarkably better that sEM. By the annealing procedure,  the clustering results should be improved much better.\n(4).  The authors claim that the proposed training scheme is  good for non-stationary, high-dimensional streaming data, but there are not much analytic results with the dimenionality and non-stationary streaming data.  The forgetting  results by controlling the annealing parameter is too rough. ", "extra_scores": {}, "extra_text": {}}, {"review_id": "a_qbJG_-RDN", "reviewer": "AnonReviewer2", "summary": "", "questions": "", "limitations": "", "rating": 5, "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 technique to training GMM using SGD instead of (s)EM. The major contributions are:\n1. a proposal for numerically stable GMM training by SGD\nthis is achieved with max-component log-likelihood approximation and an annealed process to smooth the SGD.\n2. an automatic annealing procedure that ensures SGD convergence from a wide range of initial conditions without prior knowledge of the data (e.g., no k-means initialization) which is especially beneficial for high-dimensional data\nSection 4.4 has shown that such annealing hyper-parameter can control re-learning and retention process, which is good.\nHowever, section 4.1 has shown that with or without random initialization does not impact the perf too much, then why not just use k-means initialization? is it very costly? I believe k-means initialization is also a randomized process.\n\nHowever, there are some issues:\n1. When comparing SGD vs sEM, two strong assumptions are made: \n1) SGD annealing has converged and \n2) GMM responsibilities are sharply peaked so that a single component has responsibility of around 1\nIt basically requires the data does not have a lot of noise, where each point can be assigned with to a label with a dominating probability, so what happens if the data has some noise, and how will such solution reacts for different level of noise?\nA second question is: is there any theory or bounds to support the convergence assumption? all experiments do SGD for 2 or 3 epochs, how will the loss be like after 3 epochs?\n2. No actual examples or discussions are given in terms of numerical instabilities when data dimension is high\n3. For section 4.3, the streaming scenario is not clear, is mini-batch size constant? random? if fed data of one sample, is the SGD stilling running for 2 epochs? or halting and wait?", "extra_scores": {}, "extra_text": {}}, {"review_id": "ZMWM5Diy8h2", "reviewer": "AnonReviewer4", "summary": "", "questions": "", "limitations": "", "rating": 5, "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": "The paper proposes a new approach to train GMMs using SGD under a variety of settings (streaming, concept drift, etc) addressing the issues of catastrophic forgetting, problem of parameter initialization and numerical instability.\n\nThese are all important and interesting challenges that can help advance the training of latent variable models in general. However I feel that paper is not presented with sufficient clarity to pass the ICLR bar and the exposition could be greatly improved. It is hard for me to grasp the key findings or takeaways in the paper wrt to other existing baseline methods.\n\nI also have some specific comments:\n\n- The contributions sub-section in 1.3 are vague to read, instead of simply saying - \"a novel method\", \"an automatic annealing procedure\" it would be useful to explicitly state what is novel / automatic etc - this would make it easier for the reader to understand the novelty/exact technical contributions made by the authors\n\n- In section 3.2, I am curious how the max-component approximation lower bound compares to the other lower bounds used in EM (e.g. Jensen derived lower bound or evidence lower bound in Variational Inference)? To me the presented lower bound seems like inefficient and loose (e.g. if the values are close in magnitude, then intuitively the gap between the sum and max is going to be very large). Maybe I am mistaken in understanding this bound, but it would need some more explanation and it should be contrasted with the vanilla lower bounds used in EM.\n\n- In section 3.3, the authors mention updating a subset of components to break the symmetry. I am curious how a baseline of simply perturbing the GMM parameters randomly would perform since that would also help break the symmetry?\n\n- I liked the idea of using annealing to avoid local optima and this seems to be one of the key contributions of this work in my opinion. One of my main questions here is: how did the authors measure the value addition of annealing? Did they compare the final solutions obtained by the proposed approach with a baseline (without any annealing) that uses random initialization to deal with local optima?\n\n- Notations in the paper could be improved further to make it more readable. For e.g. in Algorithm 1, the iteration steps are missing in the updates (would be good to use t=1, 2, ..T and use them in the update equations since it is not clear which iterate time steps the parameters on the RHS belong to). Also, maybe I missed this while reading but prec_clipping(..) doesn't seem to be defined near the algorithm section. \n\n- I found some vague statements in the Empirical section which raises lot of follow up questions. Would suggest making the description and analysis of the results more precise. e.g. In the beginning of section 4, authors say: \".. repeated 10 times with identical parameters..\" -> what does identical mean? what was varied?\n\n\n\n", "extra_scores": {}, "extra_text": {}}, {"review_id": "QXMhJ9K9874", "reviewer": "AnonReviewer5", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 2, "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": "A major concern about the paper is related to the unsupported claims and contribution throughout the paper. For example, the way the training copes with distribution shift or alleviate forgetting is not clear or elaborated on. Beyond the abstract and before the empirical validation no theory or justification is provided to substantiate this claim. The idea of the paper and the motivation are very interesting. The experiments look convincing. Writing and presentation are a good start point for improving the paper. \n", "extra_scores": {}, "extra_text": {}}, {"review_id": "BKluMxDCrK", "reviewer": "AnonReviewer3", "summary": "", "questions": "", "limitations": "", "rating": 5, "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": "This paper presented a stochastic gradient descent approach to learn a non-stationary high-dimensional Gaussian mixture model from online data. The authors identified 3 challenges - local optima, numerical instability, and catastrophic forgetting, and proposed to address these challenges respective with adaptive annealing, exponential-free approximation, and adaptive SGD learning rate. The proposed approach is demonstrated with several vision/non-vision tasks.\n\nOverall, I feel that the paper is slightly below the borderline. There lacks some theoretical analysis of the proposed ideas, an approach to identify the number of mixture components, and an argument as to why GMM is preferred instead of other representation learning techniques.\n\nPros:\n\n+ Interesting combination of new research trends (continual learning) and old models (GMM).\n\nCons:\n\n- Lack of an approach to identify the number of mixture components\n- Lack of theoretical justification about the max-component approximation and soft max-component approximation.\n- Lack of demonstration on why catastrophic forgetting is avoided and how nonstationary data affects this and other algorithms in experiments.\n\nSo the following improvements could improve my ratings: \n\n- An empirical analysis of how the proposed approach avoids catastrophic forgetting with nonstationary data, and a theoretical analysis/comparison between existing approaches, such as regularization, memory replay, and network morning.\n- A mathematical justification about the soft max-component approximation: how it is related the classical GMM, to what extent it is an approximation.\n- A theoretical or emprical approach to adapt the number of mixture components, for example, through introducing a prior. \n", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: This paper proposes training Gaussian mixture models using SGD, creating an algorithm appropriate for streaming data. However, we feel that the current manuscript does not sufficiently support the proposed method, and lacks insight into its workings. The reviewers believed the method lacked justification (while the authors claim to have added theoretical justification to the revised manuscript, I did not see any such new theory), and were not convinced that the method offered a significant improvement on existing methods. 
Decision: Reject