Feature learning from non-Gaussian inputs: the case of Independent Component Analysis in high dimensions
Authors: Fabiola Ricci
Lorenzo Bardone
Sebastian Goldt
Venue: ICML 2025 spotlightposter (2025)
Keywords: independent component analysis
stochastic gradient descent
FastICA
Primary area: theory
TL;DR: We establish nearly sharp thresholds for the algorithmic sample complexity of independent component analysis, an unsupervised learning algorithm that learns filters similar to deep convolutional networks.

Abstract: Deep neural networks learn structured features from complex, non-Gaussian inputs, but the mechanisms behind this process remain poorly understood. 
  Our work is motivated by the observation that the first-layer filters learnt by deep convolutional neural networks from natural images resemble those learnt by independent component analysis (ICA), a simple unsupervised method that seeks the most non-Gaussian projections of its inputs. 
  This similarity suggests that ICA provides a simple, yet principled model for studying feature learning. 
  Here, we leverage this connection to investigate the interplay between data structure and optimisation in feature learning for the most popular ICA algorithm, FastICA, and stochastic gradient descent (SGD), which is used to train deep networks. 
  We rigorously establish that FastICA requires at least $n\gtrsim d^4$ samples to recover a single non-Gaussian direction from $d$-dimensional inputs on a simple synthetic data model. We show that vanilla online SGD outperforms FastICA, and prove that the optimal sample complexity $n\gtrsim d^2$ can be reached by smoothing the loss, albeit in a data-dependent way. We finally demonstrate the existence of a search phase for FastICA on ImageNet, and discuss how the strong non-Gaussianity of said images compensates for the poor sample complexity of FastICA.

Feature learning from non-Gaussian inputs: the case of Independent Component Analysis in high dimensions
Fabiola Ricci! Lorenzo Bardone
Sebastian Goldt
Abstract
tional neural networks (CNNs), which recover Gabor fil- ters (Gabor, 1946) like the ones shown Figure I(a) in their first layer when trained on sets of natural images like ImageNet (Krizhevsky et al., 2012: Lindsey et al., 2019; Guth & Ménard, 2024). Since Gabor filters are localised in both space and frequency domains, they are effective for capturing edges and textures in natural images (Hyvärinen et al., 2009). How these filters emerge from training a deep network end-to-end on natural images with stochastic gradi- ent descent remains key question for the theory of neural networks.
Deep neural networks learn structured features from complex, Gaussian inputs, but the mech- anisms behind this process remain poorly under- stood. Our work is motivated by the observation that the first-layer filters learnt by deep convo- lutional neural networks from natural images re- semble those learnt by independent component analysis (ICA), a simple unsupervised method that seeks the most oon-Gaussian projections of its inputs. This similarity suggests that ICA pro- vides simple, yet principled model for study- ing feature learning. Here, we leverage this con- nection to investigate the interplay between data structure and optimisation in feature learning for the most popular ICA algorithm, FastICA, and stochastic gradient descent (SGD), which is used to train deep networks. We rigorously establish that FastICA requires at least Z d' samples to recover single non-Game Gan-Gassiian direction from d-dimensional inputs on simple synthetic data model. We show that vanilla online SGD outper- forms FastCCA, and prove that the optimal sample complexity n 2 d can be reached by smoothing the loss, albeit in data-dependent way. We fi- nally demonstrate the existence of a search phase for FastCCA ImageNet, and discuss how the strong non- non-Gaussianing of said images compen- sates for the poor sample complexity of FastICA.
Deep neural networks are not the only approach to learning Gabor filters directly from data. Instead, it is well-known that independent component analysis (ICA) (Comon, 1994; Bell Sejnowski, 1995; Olshausen Field, 1996; Bell Sejnowski, 1996; Hyvärinen & Oja, 2000), a simple, unsu- pervised learning algorithm, learns filters that are similar the first-layer filters of deep CNN when applied to patches of natural images, see Figure 1(b). The key idea of ICA is to identify directions of maximal non-Gaussianity in input space, suggesting that the non-Gaussian structure in data plays a crucial role in shaping the representations learnt by deep neural networks (Refinetti et al., 2023). More precisely, given set of 7ero-mean inputs D {a+1, x² ....., ...,,."), x²), ) ICA seeks the direction w² E Rd that maximises the non- Gaussianity of the projections of the inputs == W x
"" arg max ED G(w x),
=1
where the contrast function G is measure of the non- Gaussianity of the projection 8, for example its excess kurtosis G(s) = g 3. For ICA, inputs are always pre- whitened, so ICA can be seen as a refinement of stan- dard principal component analysis, which instead seeks the projection of the data that maximises its variance, i.e. G(s) = 8 On images, PCA yields a completely differ- ent of filters which spatially extended and oscillating due the approximate wranslaton-invariancee of the image patches (Hyvärinen et al., 2009), see Figure I(c).
1. Introduction
The practical success of deep neural networks is gener- ally attributed to their ability learn the most relevant input features directly from data (LeCun et al., 2015). A classic example for this feature learning are deep convolu-
'Interatioonal School of Advanced Studies. Trieste, Italy. Corre- spondence to: Fabiola Ricci <fricci @ sissa.it>, Lorenzo Bardone <Ibardone @ sissa.it>, Sebastian Goldt <sgoldt@ sissa.tt>.
A closer look at the dynamics of the CNNs trained on the standard ImageNet dataset (Deng et al., 2009) reveals further similarities between ICA and deep CNNs. In Fig- ure (dd), we plot the excess kurtosis of the first-layer pro-
Proceedings of the 4 International Conference Machine Learning. Vancouver, Canada. PMLR 267, 2025. Copyright 2025 by the author(s).

---

Feature learning from non-Gaussian inputs
(a)
(b)
(c)
30
AlexNet
excess kurtosis
ResNet18
20
De11eNet121
10
0
101
A
10a A
105
102
10A
10s
training steps
training steps
Figure 1. Independent component analysis yields similar filters as deep convolutional neural networks. (a) The first-layer convo- lutional filters extracted from an AlexNet trained on ImageNet show the hallmarks of Gabor filters (Gabor, 1946), like orientation, localisation, and serddwidh--pecfffiiig (b) Independent components found with online SGD on 100 000 64 x 64 patches extracted from ImageNet. (c) Leading principal components of the same ImageNet patches. (d) Excess kurtosis of the representations of the first convolutional layer of three deep CNNs during training on ImageNet. For each filter. we compute the dot product between its weights and patches sampled from ImageNet, then compute its excess kurtosis averaged over image patches and neurons. (e) Snapshots of a subset of convolutional filters of AlexNet at various points during training. As the excess kurtosis of representations increases, prototypes of Gabor filters emerge, which are then sharpened during the remainder of training. Full experimental details in Appendix A.I.
jections st. = - xx of ImageNet patches xz along the weights {w*} k==I....K of the K first-layer filters of three different CNNs during training, averaged over patches and neurons. If the projections sy are normally distributed, the excess kurtosis is (close to) zero, as is the case early during training when weights are close to their (uniform) initialisation. In all three networks, the -on-Gaussianity of projections increases sharply between 10 and 10 steps, before plateauing. Interestingly, the moment in which the deep CNN focuses on non- Gaussian projections is precisely the moment in which the Gabor filters form, see Figure I(c).
not the speed of convergence, but the length of the initial search phase before the algorithm recovers any trace of a signal (Bottou, 2003; Bottou & Le Cun, 2003). The focus of our work is therefore in establishing lower bounds on the sample complexity required by FastICA and SGD for escaping the search phase.
Recently, Auddy & Yuan (2024) showed that in high di- mensions, n Z d/s² samples are necessary and sufficient to ensure the existence of estimates of the non-Gaussian directions with an error up to & when given unlimited com- putational resources, while the sample complexity required for polynomial time algorithms is n 2 d2/e2. They also proposed an initialisation scheme for FastICA that reaches the optimal performance. While these results establish fun- damental statistical and computational limits of ICA. pre- cise and rigorous analysis of how the two main algorithms used in practice, FastICA and SGD, escape the search phase starting from random initialization is lacking.
Taken together, the similarity between filters obtained by deep CNNs and ICA and the dynamics of the excess kurtosis of first-layer representations of deep CNNs suggest that ICA can serve as simplified, yet principled model for studying feature learning from non- Gan-saussiin inputs.
In this work, develop quantitative theory of feature learning from non-Gaussian inputs by providing a sharp analysis of the sample complexity of two algorithms for Independent Component Analysis, namely FastCCA, the most popular ICA algorithm used in practice, and SGD, which is used to train deep neural networks.
In this paper, we derive sharp algorithmic thresholds on the sample complexity required by ICA recover non- Gaussian features of high-dimensional inputs. Our analy- sis leverages recent breakthroughs in the analysis of super- vised learning dynamics with Gaussian inputs Ben Arous et al. (2021); Damian et al. (2023); Dandi et al. (2024) for analysing the unsupervised case with non-Gaune Gaussian inputs. Our main results are as follows:
Early theoretical works independent component analysis focused mostly the speed of convergence of different algorithms for finding w" in the classical regime of fixed input dimension and large number of samples; the most popular algorithm for performing ICA, FastCCA (Hyvärinen & Oja, 2000), derives its name from the fact that it con- verges with a quadratic (Hyvarinen, 1999) rather than linear rate. However, ICA struggles in modern applications where inputs tend to be high- high-dimensional: for example, its perfor- mance is highly sensitive to its initial conditions (Zarzoso et al., 2006; Auddy & Yuan, 2024). For algorithms run- ning in high dimensions, the main bottleneck is typically
We prove that the popular FastCCA algorithm exhibits poor sample complexity, requiring a large-size batch of n 2 d to recover a hidden direction in the first step (see Section 3.3);
We prove that the sample complexity of online SGD can be reduced down to the computational threshold of n 2 d?, at the cost of fine-tuning the contrast function in a data-dependent way (see Section 3.5);
2

---

Feature learning from non-Gaussian inputs
We demonstrate that FastICA exhibits a search phase at linear sample complexity when trained on Ima- geNet patches, but recovers non-Gaussian direction at quadratic sample complexity, and we discuss how the strong non-Gaussianity of real images speeds up recovery (see Section 4).
dataset D, we initialise the weight vector randomly on the unit sphere, wo Unif(( ), and then iterate the FastICA updates for t 2 1 until convergence:
th Ep G"w--1 - x)] - ED|G" x) w
Wt /r/ 1
(4)
The main technical challenge in establishing our results dealing with the non- Gaussianity of the inputs, which adds complexity to our analysis both by introducing addi- tional terms whose statistics need to be controlled, and by impacting the intrinsic properties of the loss function.
The FastICA iteration is made of two contributions: a "gra- dient" term ()](w,,))lw=.-.. x G' (wi x); which
drives the algorithm towards a direction where the gradient vanishes, and the "regularisation" ( z) 2)] which ensures quadratic rather than the linear convergence of first-order methods. giving FastICA its name. As we will show in Proposition 1, the regularisation term is also key for efficient learning.
2. Setup
We first describe the most popular algorithm to perform ICA, called FastICA, and we introduce the classic ICA model for synthetic data which will allow to systematically test the performance of FastICA.
2.2. The ICA data model
To carefully test the sample complexity of FastICA and other ICA algorithms like SGD in a controlled setting, we first study the case of synthetic data according to a noisy version of the standard ICA data model before moving to real data in Section 4. The idea is have inputs that follow an isotropic Gaussian distribution in all but one direction that we will call the spike v E s The projection of the inputs along the spike yields a non-Gaussian distribution, so this is the direction that ICA should recover. Specifically, we assume that the data points D = {x*))² are drawn according
2.1. FastICA and contrast functions
Performing ICA generic set of inputs x E Rd drawn from data distribution P with zero mean and identity co- variance means finding unit vector w E gd-1 which is an extremum of the population loss
(2) where G : RR + R suitable "contrast function"' which measures the Gaussianttyy of the projections s w x. The two standard choices for G(s) used in practice (Pe- dregosa et al., 2011) are
((w) := Ep G(w- x),
x" =S =s =S(BBuut 24) E ER"
(5)
G(s) = -M-11 and G(s) == log cosh(as)
where z" (O, 11) is vector with white noise, while the scalar random variable vl is the latent variable for each input, and follows non-Gaussian distribution. For concrete- ness, we choose a Rademacher distribution under which Dd = =1 with equal probability. The signal-to-noise ratio 8 2 0 sets the relative strength of the non-Gaussian com- ponent along e compared to the Gaussian part =" of each input. Since whitening the inputs is standard preprocess- ing step for FastICA (Hyvärinen & Oja, 2000), we finally pre-multiply the inputs with the whitening matrix
a
for some E [1. Another classical choice is the excess kurtosis, G(s) := 54 - 3, even though it is more sensi- tive to outliers; therefore, the first two contrast functions are preferred in practice since their growth is slower than polynomial. In practice, the expectation is approximated as an average over a set of n centered and whitened inputs D (xx) 1== (that we assume to be drawn i.i.d. from IP).
As an optimisation algorithm, we focus on the classic FastICA algorithm of Hyvärinen & Oja (2000). While number of works have recently proposed alternative algo- rithms with provable guarantees, e.g. Arora et al. (2012). Anandkumar et al. (2014), or Voss et al. (2015), we focus in- stead on FastICA since it is the most popular ICA algorithm (it is the reference implementation in seikit-learn (Pedregosa et al., 1)).
8 S= 1+ 1+8- +
vv
dxd
which ensures that the covariance of the inputs simply the identity. This whitening step makes impossible to trivially detect the non- nnn-Gausiian direction by performing PCA or any other method which is based only on the information in the first and second cumulant. An alternative perspective of the ICA model is that it provides inputs where the noise is correlated, and in particular has a negative correlation with the spike to be recovered.
FastICA finds extremal points of the population loss (2) via a second-order fixed-po point iteration. Here, we simply state the FastICA algorithm for convenience; see Hyvärinen & Oja (2000) for a detailed derivation and discussion. Given a
3

---

Feature learning from non-Gaussian inputs
1.-
d²
n n=d?*t 0.2
d
0.8
Overlap
0.6
0.4
02
4.54
0.0
FastiCA steps
FastICA steps
FastICA steps
Figure 2. Performance of FastICA on the spiked cumulant model for various sample complexity regimes. We run the standard FasICAA updates (4) data drawn from noisy ICA model (5). At each step of the algorithu, we draw new batch of size df (left), = d3+ with ε = 0.2 (middle) and n = d² (right). In the latter case, we also show the corresponding loss function. We plot the overlap between the planted non- Gaussian direction and the estimate w produced by FastCCA. While recovery of the spike theoretically possible at n = dd (see Section 2.3), FastICA does not recover the spike at that sample complexity, and instcad requires n d² samples to recover the spike reliably. Theorem confirms this picture for a general data model. Full experimental details Appendix A.2 and additional cases = 0.5, 0.8 in Figure 5.
2.3. Fundamental limits of recovery
FastICA for couple of steps with n = d samples. In the plot, we consider the limit of online learning, in the sense that we draw a new training set at each step of the algorithm: in Figure 7, we show that the picture does not change when using the same training set repeatedly; see also our diseussion in Section 3.3.3.
We now briefly recall the fundamental limits of recovering the spike in the ICA data model (5), which will serve as benchmarks to evaluate the performance of specific algo- -hheorttitt
rithms like FastICA or SGD. From an information thcoretic point of view, i.e. when given unbounded computational resources, the presence of the Gaussian direction v can already be detected using a linear number of samples n 2 d. However, detecting the spike at this sample complexity is only possible through exponentially complex algorithms that are based on exhaustive search. Efficient algorithms, that can run in polynomial time in d, require at least a quadratic number of samples n = wdd).
We can see immediately that at quadratic sample complexity, FastICA does not recover the spike at all. Instead, the overlap is stuck at a value of roughly 1/vd,, which is the typical overlap of a randomly direction in IR with the spike. This small value for the overlap suggests that the algorithm is stuck in the search phase, and does not recover the spike at quadratic sample complexity. With instead, FastCCA recovers the spike within a few steps. In the intermediate regime of n = .B++ samples, & = 0.2 the algorithm recovers the spike more slowly, and there are large fluctuations.
d' samples
For details on the algorithmic result, Szekely et al. (2024); Dudeja & Hsu (2024) and Auddy & Yuan (2024) for a thorough analysis using low-degres methods, and Di- akonikolas et al. (2017) for Statistical Query (SQ) analysis, which reaches equivalent predictions using the low-degree to SQ equivalence (Brennan et al., 2021). Auddy & Yuan (2024) establishes the coarmattionthhoretie limits of recon- structing the spike using min-max type analysis. A similar statistical-to computational gap appears in the related prob- lem of tensor PCA (Richard & Montanari, 2014).
These experiments suggest that FastICA has a poor sample complexity that is quartic in the input dimension for recov- ering the spike, with intermediate regime with (0, 1). Our theoretical analysis of FastCCA will now confirm that these are indeed the fundamental limits of FastCCA.
d++
2.4. Experiment: FastICA requires a lot of data to recover the non-Gaussian direction
3. Theoretical Results
Before diving into our analysis of FastICA and stochastic gradient descent (SGD) for performing ICA, we highlight how to address the key technical difficulty in analysing learning from non-Gaussian inputs.
We are now in position to evaluate the FastCCA algorithm in the controlled setting of the noisy ICA model (5), where the goal is to recover the planted non-Gaussian direction (or "feature") v. Specifically, we can compare the sample complexity of FastICA, i.e. the number of samples FastCCA requires to recover v, to the fundamental limits of efficient or polynomial-time recovery, which is d²
3.1. ICA, likelihood ratios, and information exponents
The key difficulty in analysing the ICA loss (2) or the FastICA update (4) is that it entails an average over an explicitly non non-Gaussian distribution of the inputs. This is in
In Figure 2, we plot the absolute value of the overlap be- tween the spike and the estimate w obtained by running
4

---

Feature learning from non-Gaussian inputs
contrast the vast majority of works on (supervised) learn- ing algorithms, which usually consider a Gaussian input dis- tribution, an input distribution that is effectively Gaussian via the Gaussian Equivalence Theorem, or inputs that are distributed uniformly on the unit hypercube (as discussed at the end of the section).
The contrast function G and the data distri- bution IP are such that the loss (2) is a function of the overlap a v, for a given spike with unit norm. Moreover, we assume that f(a)
Assumption 1.
d2 (a) E LL2 (R. dPo
This assumption is trivially satisfied for the noisy ICA model (5) and the standard contrast functions in Equa- tion (3), but it allows us to prove our results on the sample complexity of FastICA (Theorem 2) and smoothed SGD (Theorem 4) for arbitrary sub-Gaussian latent variables v (see Definition in Appendix B). In principle, these hy- potheses exclude the case of fat-tailed distributions on the latent variable, i.e. those that do not belong to L (IR, Po). However, the results of simulations with a Laplace prior on V (see Figure 9) reveal the same picture as in the sub- Gaussian case, suggesting that this assumption could be further relaxed.
The key idea of the theoretical analysis of FastICA and SGD is to rewrite the average over the non-Gaussian in- puts average over the standard normal distribution Po = (0, Ia), as
(ww)= == Er (w x = Er, [G x) x)]
(6)
by introducing the likelihood ratio
dP (s):
((s)
(7)
akin to how one would proceed in analysis with the second moment method in hypothesis testing (Kunisky et al., 2019). Intuitively, the likelihood ratio (or more precisely, its norm) measures how different the distribution of the projection s = w. is from a standard Gaussian distribution. The likelihood ratio will thus be the key object in our study to determine the number of samples that a given algorithm like FastICA requires to find projection that is distinctly non-Gaussian and yields large likelihood ratio.
Assumption 2 (Contrast function). We assume G to be even. For our analysis of FastICA, we require G E (RR) with bounded derivatives. For smoothed SGD in Sec- tion 3.5, we additionally require Σ Exzo(cE)) < and that there exist some C1, C2 > 0 such that, for any E E. [G'(s) S C1 + s2)Ca holds.
Under the noisy ICA model (5), the likelihood ratio (7) depends only on the projection of the inputs along the non- Gaussian direction v to be learnt. We can expand the popu- lation loss in terms of single scalar variable, the overlap := w- e between weight w and spike v:
3.3. Sample complexity of FastICA
Our goal is to rigorously quantify the sample complexity of FastICA on the noisy ICA model (5). FastICA proceeds by performing few iterations of the update (4) with the full data set. This setting is precisely the large-batch (Ba et al., 2022; Damian et al., 2022) or "giant step" (Dandi et al., 2024) regime that has attracted a lot of theoretical interest recently. Following the approach of these recent works, the key idea of our analysis is to study whether a single step of FastCCA with a large number of samples x d yields an estimate wi that has an overlap with the spike e that does not vanish with the input dimension.
E k!
((w)= = BF. GGw x)E(V x)] = k20
a
(8)
where the last equality has been obtained by expanding the contrast function G and the likelihood ratio l in series of Hermite polynomials and by using the orthogonality of Her- mite polynomials (Lemma 11 in Appendix B.3). A key quantity for the analysis of ICA is the order k* of the first term in the expansion (8) that has non-zero coefficient. This quantity was introduced as the "information exponent" by Ben Arous et al. (2021) in the context of supervised learning, when they showed that it governs the sample com- plexity of online SGD in single-index models, like the noisy ICA model. Related quantities like the leap index of Abbé et al. (2023); Dandi et al. (2024) and the generative expo- nent of Damian et al. (2024) govern the sample complexity for more complex data models with several directions to be learnt. Here, we note that due to the whitening of the inputs, the population loss (2) for the noisy ICA model model has information exponent A* 4, as we show in Appendix C.
1. THE IMPORTANCE OF REGULARISATION
It is instructive to first consider briefly the ideal case of infi- nite data. The following Proposition shows that FastCC is able to fully recover the spike v as d + oo if we re- place the empirical average Ep with the population average Ep. More interestingly. we show that it is the regularisation term, originally introduced to speed up convergence (see Section 2. 2), that gives the essential contribution to escape mediocrity at initialisation.
Proposition 1 (Infinite batch size). Assume that the inputs are distributed according to the noisy ICA model (5) and that of + 0- this holds for the standard contrast functions in Equation (3), Set a1 == W1 V, where we is the updated weight vector given by the first iteration of FastCCA.
3.2. Main assumptions
For the following theoretical analysis, we assume that:
5

---

Feature learning from non-Gaussian inputs
Then, in the infinite batch-size limit, we have that
overlap in that regime. For quadratic sample complexity, the gradient in the FastICA update does not concentrate, so the algorithm does not recover the signal and remains stuck in the search phase.
of =
On the other hand, if the regularisation term is missing and only the "gradient step" is taken, i.e. the iteration reads We = Flxx''uu .)), then we obtain of = O(1/))
Key elements of the proof. The overall ideas of the proof are similar to those employed in the analysis of the teacher- student setup of Theorems and 2 from (Dandi et al., 2024), with the following main differences. First of all, in our unsu- pervised case, instead of isotropic inputs, the data points are in general non- non-Gaussian distributed; as already observed in Section 2.3, this impacts the information exponent k- of the population loss. Second, unlike for SGD, the iteration of FastICA presents both a gradient term and a regularisation term which is data-dependent (the regularisation "constant" depends on G" (w z),, and its analysis is hence just as com- plex as that of the gradient term). The general strategy is computing the expectations for both the signal Wi v and the noise || ||. and then identifying how much data n al- lows for the concentration of the two quantities. To be able to compute the expectations, we exploit the orthogonality properties of Hermite polynomials proved in Corollary 14, Once concentration of both signal and noise is established, the positive result corresponds to the signal dominating the noise, and vice versa for the intermediate regime in the mid- dle. Conversely, if the signal does not concentrate, is not sufficiently large for the overlap escape from its scaling at initialisation and first upper bound holds.
3.3.2. FINITE SAMPLES: HARDNESS OF LEARNING
We now provide almost sharp bounds on the sample com- plexity of FasICCA for learning the Gan-Saasiinn feature v. For the analysis, we follow the recent "one-step" analy- ses (Ba et al., 2022; Damian et al., 2022; Dandi et al., 2024) and study the overlap of the spike with the weight obtained from performing a single step of FastICA with a large batch of n = dh samples. More precisely, we provide a lower bound (positive results) for the amount of signal that can be weakly recovered given number of samples on the order of d² where k* is the information exponent of the loss, and precisely& k* = 4 for the noisy ICA model. In the data-scarce regime where n = o(d ), we provide upper bounds (nega- tive results) on the overlap a = w . v that can be achieved in a single step of FasuCA. We distinguish two cases: either there is no improvement with respect to the random initiali- sation, or the upper bound on the amount of signal learned is cimeenioonality-deeeddee
meaning that n O(d' -1) are still not sufficient to exit from the search phase in the -igh-dimensional limit. Yet this regime smoothly bridges the situation of total ignorance, where nothing is learnt at all, and recovery of the spike in one step, thanks to a continuous dependence on the parameter 5.
Theorem 2. (Finite batch size) Let k- be the information exponent of the population loss (2). Consider a number of
3.3.3. DISCUSSION
samples n ((dk'- for å E [0, 2). Set ai == W1 c, where wi is the updated weight vector given by the first iteration of FastICA. Then,
The almost sharp analysis of the high-dimensional asymp- totics of FastCAA ih 2 is our first main result, and explains theoretically why FastICA tends to struggle in high dimensions (Zarzoso et al., 2006; Auddy & Yuan, 2024). While FastICA is often applied after reducing the dimension of the inputs using PCA, it should be noted that the scaling of n x d on the noisy ICA model poor even in moderate to small dimensions. One non-standard way to escape this scaling is to use the ad-hoc initialization scheme proposed recently by Auddy & Yuan (2024), which allows to find an initialisation that gives ICA a warm start with non-trivial overlap given only n 2 d² samples. However, there remains large gap between the performance of FastICA as it is com- monly used, and the 2 d² or, more generally speaking, n 2 d2 /2 bounds we have from low-degree methods, see (Damian et al., 2024). One may therefore ask whether other algorithms, and in particular the simple stochastic gradi- ent descent, may be more efficient at performing ICA, not least since neural networks trained with stochastic gradient descent seem to learn ICA-like filters in their first layer rapidly.
E (1,2) 8 E (0,1]
+ of = 0 (t)
of = 0
0=0 0
of 21 1 -().
We prove Theorem 2 in Appendix D. Even if Theorem 2 is not specific to the noisy ICA model, and instead applies to any input model and contrast function satisfying Assump- tion and 2, it is instructive to look at the special case of Rademacher prior on the latent variable and standard con- trast function. Equation (3). In this setting, we have proved in Appendix C. in the appendix that k* 4. Therefore, we find that n O(d*) is sufficient to perfectly recover the spike with single step. up to contributions which vanish when d goes to infinity, in line with what we see in the ex- periments of Figure 2. Moreover, if d S n & d², it turns out that the overlap a1 = w1= U is bounded by a dimension- dependent constant, which explains the slow increase of the
6

---

Feature learning from non-Gaussian inputs
3.4. Vanilla SGD
Biroli et al. (2020) introduced the idea of smoothing in the context of tensor PCA. Damian al. (2023) showed that smoothing the loss for a single neuron ow x) with activation function o : R + R trained on supervised task with Gaussian inputs and labels provided by a teacher neuron a(v x) reduces the sample complexity of weak recovery for online SGD from d -1 (Ben Arous et al., 2021) down to the optimal d /2, Their analysis crucially relies the fact that teacher and student have the same activation function For ICA, one can roughly think of the contrast function as the student activation function, which will in general have a different information exponent from the loss function, which depends also on the likelihood ratio. We therefore generalise the analysis of Damian et al. (2023) to this mismatched case and show that the dynamies of smoothed SGD on ICA is governed by the interplay of the information exponents of the population loss and the contrast function respectively, called kij and kj, and that the optimal sample complexity is reached if and only if kj = k2- This fact follows immediately from the heuristic analysis performed in Section 3.5.1, where it is shown that for the optimal X = dd'4, the ODE describing the behavior of the overlap a at the beginning of learning is
In order take a closer look at feature extraction in the context of deep neural networks, we analyse SGD alternative algorithm to perform ICA. Consider set of centered and whitened data points DD = {a") J/=-- We sample a new data point at each step. For a suitable learning rate d > 0, each iteration of the spherical online SGD, for wo ~ Unif((d-1) 1) is defined as
W = - + 8 V vpn C(u 11) WEE
t 121 1,
w / well,
The spherical gradient sph for a function f is given by V su, := (Ia ww V (ur.). In the case of the noisy ICA model model (5), we directly apply the analysis performed by Ben Arous et al. (2021), which guarantees that if n = log" d), or more in general n = 21d log" d), the spike is "strongly" recovered, meaning that wn v - in probability for n going to infin- ity. Moreover, n- = Odd)), the spike cannot be recovered, i.c. supiga/o 1u wel + 0 in probability.
Hence, vanilla online SGD is already faster that FastICA in the large-batch setting. However, given the information- theoretic results of Section 2.3, online SGD does not achieve the optimal performance and there remains a statistical-to- computational gap.
a(t) di k5/2
d'(()
(9)
On the noisy ICA model (5), for which kj = 4, we thus find the following picture. For the standard contrast functions (Equation (3)), ky * k5 = 2, and then smoothed online SGD does not improve the sample complexity of vanilla online SGD: the spike still recovered in a cubic time. In contrast, if k5 = 4, i.e. the information exponents of the loss and the contrast function are equal, smoothed online SGD is able to recover the spike in a quadratic time, which means that smoothing the landscape reduces the sample complexity of online SGD down to the optimal threshold suggested by low-degree methods. Even if it would be clearly necessary a fine-tuning of the contrast function in a data-dependent way which cannot be performed in practice, we observe that in our case the fourth-order Hermite poly- nomial G(s) ha(s), which for whitening data is nothing but the excess kurtosis G(s) = s² - 3, is indeed an optimal contrast functions. Remarkably, the kurtosis is the most classical measure of non- non-Gaussinnity that one can find in the literature (see e.g. Hyvärinen & Oja (2000)) and has been employed study recent developments of FastCCA in (Auddy & Yuan, 2024).
3.5. Closing the statisticll-o-ccompuuattioonnl gap: smoothing the landscape
We now show that the statistical-to- co-computational gap of SGD on noisy ICA can be closed by running online spherical SGD on a smoothed loss following the approach of Biroli et al. (2020); Damian et al. (2023):
Definition
(Biroli el al., 2020; Damian et al., 2023) For x= 2 0 and x. E Rd,, the smoothed loss function is LA(w, (w, 2) == Lx[G(w - 2)], where the smoothing opera- tor LA reads
[o(
Au w
]
[x[G(w-x) == Expppi
Az|l
if pm is the uniform distribution over Sd-1 conditioned on being orthogonal to w.
The key intuition behind the smoothing operator is that, thanks to large (than scales with the input dimension), it allows to evaluate the loss function in regions that are far from the iterate weight w,, collecting non-local signal that alleviates the flatness of the saddle around = 0 of L, re- ducing the length of the search phase; we illustrate the effect of smoothing on the loss in Figure 6. For implementation purposes, we give an explicit formula for the smoothed gra- dient of said loss in Appendix E. Note that the "large steps" " allowed by a large X help in escaping the search phase.
1. HEURISTIC DERIVATION
We first heuristically derive the sample complexities of The- orem 4, following Damian et al. (2023). We will use Corol- lary 13 to deal with expectations of multiple products of Hermite polynomials appearing due to non-Gaussian inputs.
7

---

Feature learning from non-Gaussian inputs
1,0
Smeothed online SGD
Smoothed online SGD
0.8
Vanilla online SGD
Vanilla online SGD
Overiap
Smoothed online SGD Vanilla online SGD
0.2
0.0
n/d
n/d
n/d3
Figure Performances of smoothed online SGD and vanilla SGD on the noisy ICA model, for G(s) ha(s) and G(s) exp(-s2/2) respectively. Left: In a linear time, nor smoothed and vanilla SGD recover the spike. Middle: In the quadratic regime, smoothed SGD recovers, although vanilla SGD does not. Right: In the cubic regime, both smoothed vanilla SGD recover the spike. Experimental details in Appendix A.3.
A crucial quantity is the signal-to-noiss SNR since, as shown in Damian al. (2023), at the beginning of leaning the over- lap == w v is given by the ODE a' a1) = /NR/o(t) with
Since 2' ((I/Vd), can avoid the dependence the replicas by simply choosing À < d//4, This choice corresponds to Equation (15) in Appendix E. Hence, = (() We can conclude that the noise scales as
SNR- = Egg- u/²
E 9
Ellg |3] = (dX 2
where is the online spherical gradient g == V sph (x, w). We compute how the SNR scales with the dimension of the inputs d.
By combining signal and noise, we have that
a² 12(kg-2)
SNR
A d (1I-1) d
d -1
Signal Recall the scaling computed by Damian et al. (2023) for the signal in the case of even contrast function:
This implies that, in the case of the optimal X d the ODE describing the behavior of the overlap a for small 2 0 is nothing but Equation (9).
Elg w) = O(a d (4g (d (41- /2²)
Therefore, the number of data points required to reach an order-one overlap is d /2, i.c. the optimal one accord- ing to the fundamental limits of recovery, and only kj = kū and than we have heuristically derived the results of Theorem 4. The following theorem, which we prove in Appendix makes rigorous the heuristics. It holds for any input model satisfying Assumption Assumption and the following Assumption 3:
when a 1/2 Note that it depends the information exponent Kj of the loss.
Noise
We look at the scaling of the noise B[[bg[]?. By definition, O (>~tx (G (w (w-<x)) xr))) Since with high probability | = O(d), we study the second moment of = Cx (w z)). By exploiting the definition of the smoothing operator and writing the expectation in terms of the likelihood ratio, we get that, for two replicas 2, ~ pas
Po) Assumption 3. We assume & E 1, (Rd, Po) and that the loss (2) is strictly monotonically increasing in a = w v.
+ A=
w Az' 20 E(v x))
I]
Theorem 4 (Escaping mediocrity). Define kj and k2 as the information exponents of the population loss (2) and the con- trast function, respectively. Assume that Wy ~ Unif(SA-1) and that the initial overlap ag == wb v such that an 2 1/vd. Consider X E [1, Then, there exist n7 R and E N. satisfying
E
R
E
G"
G
V1
V1+
2.2'
E
E
G²
+ A=
x
w
+ d='
vI
x
x
By expanding G' and L, we get a Gaussian expectation of triple product of Hermite polynomials, whose scaling can be computed thanks Lemma 12, recalling that G is even and then e = 0. In particular, thanks to Corollary 13, we have
n - polylog(d)).
17 O 2/2 polylog(d))
such that On 2 1 - 1-4 with high probability. Here,
On == Wn v and Wn the estimator given by smoothed online SGD with Th = n) and Xn. X.
EP, (G'(w) r)G'(wz x)((v= x)] = @((w) )
8

---

Feature learning from non-Gaussian inputs
(a)
n = 2d
(b)
11 = 2d
(c)
2.5
0.00
0.00
eto (c(00)
ctt) (c(O))
2.0
-0.02
-0.02
2d
dimension d
dimension d
1.5
3=
1024
36
100
4d
-0.04
256
2304
-0.04
64
144
1.0
5=
0
0
2
20
60
steps
steps
width of ImageNet patches
Figure 4. Performance of FastICA on ImageNet. We plot the difference between the test loss [tt at step t, and the test loss (C())) initialisation averaged over twenty initialisations. At linear sample complexity (a), the test loss remains close to its initial values, but it clearly decreases at quadratic sample complexity (b). (c) Ratio of the fourth moment of projections s, of ImageNet patches of the given width along Gabor filter obtained via FastCCA, and of the fourth moment of projections so along a random direction. Different colours refer to the size of the training set from which the Gabor filters were obtained. Full experimental details in Appendix A.4.
4. The search phase of CastICA on images
Indeed, we found strong finite-size effects when attempting to measure the effective signal-to-noise ratio of the non- Gaussian directions in the images. In the spiked cumulant model, an empirical measure of signal-to-noise ratio can be obtained by dividing the fourth moment of the projection of inputs along the spike, 8. = v-z, with the fourth moment of inputs projected along a random direction so = Wi- x. As a surrogate for the "spike" in real images, we ran FasuCA to obtain clean Gabor filters from n 2, 4, 5d ImageNet patches (as detailed in Appendix A.2) and computed the ratio of the fourth moment of projections s. along these filters with projections 8o along random directions. For the spiked cumulant with standard Laplace prior, this ratio is equal to 2. In Figure 4(C), we see that this ratio for images tends towards an larger value, but we note that the empir- ical estimate only converges slowly with width, hinting at important finite-size fluctuations.
We finally investigate whether FastICA also exhibits an extended search phase on real data. To that end, we ran FastICA on patches sampled at random from the ImageNet data set (Deng et al., 2009). In Figure 4, we show the difference between the test loss [(t) at step I of FastICA and the loss (C())) averaged over twenty randomly chosen initial conditions. We find that at linear sample complexity when n = 2d, Figure 4(a), FastICA is clearly stuck in a search phase: the loss does not deviate significantly from its value at initialisation. Meanwhile in the quadratic regime at 2dd, Figure 4(b). FastICA does recover non non-trivial Gaussian direction from the start, faster than our theory would suggest.
Where does this non- non-Sausiian signal come from that FastICA picks up in the images? A natural idea would be to suggest that comes from the third-order cumulant, which non-zero for real images contrast to the spiked cumulant with Rademacher latent v. The corresponding information exponent K* would then suggest that dh- 11 is the transitionary regime, where recovery is possible in finite dimensions. However, since the contrast functions used in FastCA are symmetric, FastICA is blind to information carried by odd cumulants. We demonstrate this explicitly by constructing latent distribution which has mean zero, unit variance, a non-trivial third-order cumulant, but zero fourth-order cumulant, and find that FastCAA does not pick up any signal, see Figure& 8. Another potential explanation is that the latent variables of images are super-exponential, and follow an approximate Laplace distribution (Wainwright & Simoncelli, 1999). However, we verified experimentally that when the latent variable w² in the noisy ICA model has a Laplace distribution, Theorem 2 still accurately describes the sample complexity required for recovery, see Figure 9. However, the longer tail induces stronger finite-size fluctua- tions, which might contribute to the success of FastCCA at quadratic sample complexities.
5. Concluding perspectives
Despite their simplicity, Independent Component Analysis yields similar filters as the early layers of deep convolutional neural networks, and therefore offers an interesting model to study feature learning from non-Gaue non-Gaussian inputs. Here, considered the problem of recovering single feature en- coded in the non- Gan-Gassii input fluctuations and established almost sharp sample complexity thresholds for FastCCA and SGD. Our analysis revealed the poor sample complexity of FastICA, which might explain some of its problems in high-dimensional settings, while our experiments suggest that images have strong enough non- non-Gaussian features to compensate for this Going forward, it will be intriguing to extend our analysis to models with several non- Gaussas di- rections, along the lines of Dandi et al. (2024) and Bardone & Goldt (2024). Finally, developing model for synthetic data from which localised, oriented filters can be learnt looms as an intriguing challenge (Ingrosso & Goldt, 2022; Lufkin al., 2024).
9

---

Feature learning from non-Gaussian inputs
Impact statement
K., Weller, A., Oliver, N., Scarlett, J., and Berkenkamp, F. (eds.), Proceedings of the 4 International Confer- ence Machine Learning, volume 235 of Proceedings of Machine Learning Research, pp. 3024-3045. PMLR, 2024.
This paper presents work whose goal is to advance the field of Machine Learning. There are many potential societal consequences of our work, none which we feel must be specifically highlighted here.
Bell, A. and Sejnowski, T. J. Edges are the 'independent components'of natural scenes. In Mozer, M., Jordan, M., and Petsche, T. (eds.), Advances in Neural Information Processing Systems, volume 9. MIT Press, 1996.
Acknowledgements
We thank Mitya Chklovskii for valuable discussions. SG acknowledges funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme, Grant agreement ID 101 166056, and funding from the European NexioneeeGGGnrtttiinnEUU in the framework of the PRIN Project SELF-MADE (code 2022E3WYTY - CUP GS3D23000780001,, SG and FR ac- knowledge funding from Next Generation EU, in the context of the National Recovery and Resilience Plan, Investment PE1 Project FAIR "Future Artificial Intelligence Research" (CUP G53C22000400066.
Bell, A. J. and Sejnowski, T. J. An information- maximization approach to blind separation and blind deconvolution. Neural computation, 7(6): 129-1159, 1995.
Ben Arous, G., Gheissari, R., and Jagannath, A. Online stochastic gradient descent non-convex losses from high-dimensional inference. J. Mach. Learn. Res., 22(1), 2021.
Biroli, G., Cammarota, C., and Ricci- Tersenihrs,,,h,, F. How to iron out rough landscapes and get optimal performances: averaged gradient descent and its application to tensor pea. Journal of Physics A: Mathematical and Theoretical, 53(17):174003, 2020.
References
Abbé, E., Adserå, E. B., and Misiakiewicz, T. Sgd learning on neural networks: leap complexity and saddle-to-saddle dynamics. In Neu, G. and Rosasco, L. (eds.). Proceedings of Thirty Sixth Conference on Learning Theory, volume 195 of Proceedings of Machine Learning Research, pp. 2552-2623, PMLR, 2023. URL https: http:t:://rrpeeeeeeeeee http//peeednngssseeeeeeeeeeessee ngs. mlr malr.press) w195/abbe23a.hi .tmml.
Bottou, L. Stochastic learning. In Summer School on Ma- chine Learning, pp. 146-168. Springer, 2003.
Bottou, L. and Le Cun, Y. Large scale online learning. Advances in neural information processing systems, 16, 2003.
Anandkumar, A., Ge, R., Hsu, D., Kakade, S. M., and Telgarsky, M. Tensor decompositions for learning latent variable models. Journal of machine learning research, 15:2773-2832, 2014.
Brennan, M. S., Bresler, G., Hopkins, S., Li, J., and Schramm, T. Statistical query algorithms and low de- gree tests are almost equivalent. In Belkin, M. and Kpotufe, S. (eds.), Proceedings of Thirty Fourth Con- ference on Learning Theory, volume 134 of Proceedings of Machine Learning Research, pp. 774-774. PMLR, 2021. URL URL https:/ pttttp:s//prcceeiinge http///eeeeinggseepeeepee/e mlr mlr.pp pprsss/ v1 34 1/34/bremnan21e .tml.
Arora, S., Ge, R., Moitra, A., and Sachdeva, S. Provable ica with unknown gaussian noise, with implications for gaussian mixtures and autoencoders. Advances in Neural Information Processing Systems, 25, 2012.
Auddy, A. and Yuan, M. Large dimensional independent component analysis: Statistical optimality and compu- tational tractability. to appear in Annals of Statistics, 2024.
Buet-Golfouse.
A multinomial theorem for her- mite polynomials and financial applications. Applied Mathematics A Journal of Chinese Universities Se- ries B, 06: 00:7017 2015. URL https://api. hwaemm/ccccccoaa..ooo////ooppppp
Ba, J., Erdogdu, M. A., Suzuki, T., Wang, Z.., Wu, D., and Yang, G. High-dimensional asymptotics of feature leam- ing: How one gradient step improves the representation. Advances in Neural Information Processing Systems, 35: 37932-37946, 2022.
Comon, P. Independent component analysis, a new concept?
Signal processing, 36(3):287-314, 1994.
Damian, A., Lee, J., and Soltanolkotabi, M. Neural net- works can learn representations with gradient descent. In Conference on Learning Theory, pp. 5413-5452. PMLR, 2022.
Bardone, L. and Goldt, S. Sliding down the stairs: How correlated latent variables accelerate learning with neu- ral networks. Salakhutdinov, R., Kolter, Z., Heller,
10

---

Feature learning from non-Gaussian inputs
Damian, A., Nichani, E., Ge, R., and Lee, J. D. Smoothing the landscape boosts the signal for SGD: Optimal sample complexity for learning single index models. In Thirty- seventh Conference on Neural Information Processing Systems, 2023. URL. https: http:///ppeeeeeieeeeeeeeeee net./
Hyvarinen, A. Fast and robust fixed-point algorithms for independent component analysis. IEEE transactions Neural Networks, 10(3):626-634, 1999.
Hyvärinen, A. and Oja, E. Independent component analysis: algorithms and applications. Neural networks, 13(4-5); 41 1 1300, 2000.
Damian, A., Pillaud-Vivien, L., Lee, D., and Bruna, J. The computational complexity of learning gaussian single- 2024.
Hyvärinen, A., Hurri, J., and Hoyer, P. O. Natural image statistics: A probabilistic approach to early computa- tional vision., volume 39. Springer Science & Business Media, 2009.
index models. arXiv preprint arXiv:2403.055
Dandi, Y., Krzakala, F., Loureiro, B., Pesce, L., and Stephan, L. How two-layer neural networks learn, one (giant) step at time. Journal of Machine Learning Research, 25 (349):1-6) 65, 2024.
Hyvärinen, A., Karhunen, J., and Oja, E. In- dependent Component Analysis. John Wiley Sons, Lid, 2001. ISBN 9780471221319, doi: https///oo.ogg//0000000001111113
Deng. J., Dong, W., Socher, R., Li, L.-J,, Li, K., and Fei-Fei, L. Imagenet: A large-scale hierarchical image database. In 2009 IEEE conference on computer vision and pattern recognition, pp. 248-255. Ieee, 2009.
URL
https://onnlieellbrrryyyyy abs/10.1002/0471220 htt.--00000012222113aaaaaaeaaaaa fnatter_intttttsss
com/doi/
Ingrosso, A. and Goldt, S. Data-driven emergence of con-
Diakonikolas, I., Kane, D. M., and Stewart, A. Sta- tistical query lower bounds for robust estimation of high-dimensional gaussians and gaussian mixtures. In 2017 IEEE 58th Annual Symposium Foundations of Computer Science (FOCS), pp. 73 84, 2017. doi: 10.1 109/FOCS.2017...
volutional structure in neural networks. Proceedings of the National Academy of Sciences, 119(40))e221185411 2022.
Krizhevsky, A., Sutskever, I., and Hinton, G. E. Imagenet classification with deep convolutional neural networks. Advances in neural information processing systems, 25, 2012.
Dudeja, R. and Hsu, D. Statistical-come computaiool trade-offs in tensor PCA and related problems via communication complexity. The Annals of Statistics, 52(1):131 - 156, 2024. doi: 10.1214/23-AS22331 URL https: /doi. org/ 10. 1214 122122222222223...
Kunisky, D., Wein, A. S., and Bandeira, A. S. Notes on computational hardness of hypothesis testing: Predictions using the low-degree likelihood ratio. ISAAC Congress (International Society for Analysis, its Applications and Computation), pp. 1-50. Springer, 2019.
Gabor, Theory of communication, part The analysis of information. Journal of the Institution of Electrical Engineers-part III: radio and communication engineer- ing, 93(22)):22 1946.
LeCun, Y., Bengio, Y., and Hinton, G. Deep learning. Na- lure, 521(7553) 2015.
Guth, F. and Ménard, B. On the universality of neural encod- ings in cnns. In ICLR 2024 Workshop on Representational Alignment, 2024.
Lindsey, J., Ocko, A., Ganguli, S., and Deny, S. The effects of neural resource constraints early visual rep- resentations. In International Conference Learning Representations, 2019. URL https::/p/// /pperreveee/eeeeeeeeeeeee pet/forunte (atoormmmmd-Sxx333355555
He, K., Zhang, X., Ren, S., and Sun, J. Deep residual learning for image recognition. In 2016 IEEE Conference on Computer Vision and Pattrnn Recognition, CVPR 2016, Las Vegas, NV, USA, June 27-30, 2016, pp. 770-778. IEEE Computer Society, 2016. doi: 10.1109/CVPR.2016 90.
Lufkin, L., Saxe, A., and Grant, E. Nonlinear dynamics of localization in neural receptive fields. In Thirty-eighth Conference on Neural Information Processing Systems. NeurlPS, 2024.
Huang, G., Liu, Z., van der Maaten, L., and Weinberger, K. Q. Densely connected convolutional networks. In 2017 IEEE Conference on Computer Vision and Pattern Recognition, CVPR 2017, Honolulu, HI, USA, July 21-26, 2017, pp. 2261-2269. IEEE Computer Society, 2017. doi: 10.1 10.1109/CVPPR.2077..44.
McCullagh, P. Tensor methods in statistics. Courier Dover Publications, 2018.
Olshausen, B. A. and Field, D. Emergence of simple-cell receptive field properties by learning a sparse code for natural images. Nature, 38 (6583) 1996.

---

Feature learning from non-Gaussian inputs
Pedregosa, F., Varoquaux, G., Gramfort, A., Michel, V., Thirion, B., Grisel, O., Blondel, M., Prettenhofer, P. Weiss, R., Dubourg, V., Vanderplas, J., Passos, A., Cour- napeau, D., Brucher, M., Perrot, M.. and Duchesnay, E. Scikit-learn: Machine learning in Python. Journal of Machine Learning Research, 12:2825-2830, 2011.
Refinetti, M., Ingrosso, A., and Goldt, S. Neural networks trained with sgd learn distributions of increasing complex- ity. In International Conference on Machine Learning, pp. 11843-28863 PMLR, 2023.
Richard, E. and Montanari, A. A statistical model for tensor pea. Advances in neural information processing systems, 27, 2014.
Szekely, E., Bardone, L., Gerace, F., and Goldt, S. Learn- ing from higher-order correlations, efficiently: hypothesis tests, random features, and neural networks. In The Thirty- eighth Annual Conference on Neural Information Pro- cessing Systems, 2024. URL https:: / //ppr/eeee..e.e.eee net fett/fozee
Vershynin, R. High-dimensional probability: An introduc- tion with applications in data science, volume 47, Cam- bridge University Press, 2018
Voss, J. R., Belkin, M., and Rademacher, L. A pseudo- euclidean iteration for optimal recovery in noisy ica. Advances in neural information processing systems, 28, 2015.
Wainwright, M. J. and Simoncelli, E. Scale mixtures of gaussians and the statistics of natural images. Advances in neural information processing systems, 12, 1999.
Zarzoso, V., Comon, P., and Kallel, M. How fast is fastica? In 2006 14th European Signal Processing Conference, pp. 1-5. IEEE, 2006.
12

---

Feature learning from non-Gaussian inputs
A. Experimental details
In this appendix, we collect detailed information on how we ran the various experiments of this paper.
A.1. Figure
(a) We plot the filters of the first convolutional layer of an AlexNet trained on ImageNet using mini-batch size of 128, momentum of 0.9, weight decay of 10 10-4, cosine learning rate schedule with initial learning rate of 10-2, We also trained DenseNet121 (Huang al., 2017) and a ResNet18 (He et al., 2016) using the same recipe, but with larger initial learning rate 10-1
(b) The Gabor filters were obtained with online SGD performed 100 000 patches of size 64 64 extracted from ImageNet. We have decided to perform ICA with SGD to use the typical algorithm employed to train deep neural networks. We perform PCA before running ICA by keeping 600 principal components. The contrast function was G(s) = - exp((=222)
(c) We applied principal component analysis to the same ImageNet patches to obtain the leading principal components, and plot the eight leading ones from left to right, and top bottom.
(d) For all three networks trained on ImageNet, we took 50 logarithmically spaced snapshots of the weights during training. For each snapshot, we extracted the weights for each first-layer convolution and computed the dot-products sy. w o0!" between the weight of the kth convolution for randomly sampled patches from ImageNet wt, We normalised the & for each neuron (fixed k) to have variance I, and computed the excess kurtosis of the resulting s/. We then averaged over all the neurons to obtain the curves shown in Figure 1(b).
(c) The snapshots were taken from one of the AlexNet models.
A.2. Figure 2
Here we show the performance of FastICA on the spike cumulant model (5) in three different regimes, in the case of d = 50 and 8 = 15. The spike u has been sampled uniformly in the sphere Sd-1 The regimes correspond to the batch-sizes of n = d, n = d and n = d', respectively, for = 0.2. Each batch contains n fresh data points sampled from the distribution of the inputs IP. The contrast function G(s) = - exp(-= -2 /2). We display on the y-axis the overlap wt o for t 0, . 4 (on the z-axis), where we is the updated weight vector provided by FastICA in the first four step, plus initialisation wip gd-1 under the constraint that wi- 1/Vd. An average over 15 runs has been taken.
A.3. Figure
Here we show the velocity for recovering the spike of smoothed online SGD in three different regimes: linear, quadratic and cubic. The noisy ICA model is considered. We use d = 40, 3 15 and total number of data points of n 5d3 The parameter À is the optimal one, i.e. X d The learning rate is 17 2.9d k1/2) X225 -2), where kj 4 and kg 2, since G(s) - exp( -2 2-/2) /2). The average is made over 30 runs.
A.4. Figure
(a) and (b) We perform ICA on patches of ImageNet images with increasing dimensions, and randomly drawn. The chosen contrast function is G(s) exp(-s s² /2). The standard centering and whitening pre-processing procedure is performed before running the algorithm.
(c) We first ran FastICA on n = 2, 3, 4, 5d ImageNet image patches to obtain clean Gabor filters, which we used as surrogates of the spike v. To that end, we used the "canonical preprocessing" of Hyvärinen et al. (2009), meaning we whitened inputs and ran the FastICA algorithm only the subspace spanned by the top principal components (here, A = 50 for all widths). We then computed the fourth moment of projections s, of an independent test set of ImageNet images along these Gabor filters, and divided this by the fourth moment of projections along random directions sp-
A.5. The "3not4" prior on the latent variable, Figure 8
For this experiment, we run online FastICA as in Figure on a noisy spiked cumulant model, but with a prior over the latent variable that has mean zero, unit variance, non-trivial third-order cumulant, but zero fourth-order cumulant. More
13

---

Feature learning from non-Gaussian inputs
specifically, consider latent that take the values -1, 0 and with probabilities 1/a, 2/2 and 1/0.
B. Properties and definitions
B.1. Notation
Throughout we use the standard notations o, O, O S, w for asymptotics. We recall them here for sequences, but they generalize in the straight forward way in the case of functions.
Let (aR) KER be two real valued sequences. Then:
ag o(bb)
lim 0
k000 by
Qg E (br)
= OC > 0 R
Vk > ko
|ax S Clbsl
ag (bb))
= 3C,, C2 >
Wk > kg
C₁b² S ax I Czbk
ag e (bbb))
= 3C> 0ER
Vk ko
|axl 2 Clbbal
ak
ag E (bb))
lim
k+000
bb
Occasionally we use the shorthands an < bn for on = o(bx), an S bn for an = O(bx) and On x bn for On = (bb,).
B.2. Probability notions
We recall the properties of the space of squared integrable functions with respect to the standard normal distribution Po: Definition (Gaussian scalar product). Given f.g IRd + R. we define the 1.2-ccllr product
(f,9),, E [f(x)g(x)].
The space (RR, Po) contains all the measurable functionss f Rdd R such that ||]1| (1SSpo <
Note that this space with the product defined above is an Hilbert space, and the set of the Hermite polynomials (see next section) an orthogonal basis.
Definition 6 (Sub-Gard Sub-Gaussianity). Let Z be real valued random variable, we say that Z is nub-Gaussian if it has finite mub-Gaussian norm, defined as:
1Zlle/2 == inf {t 0 E(exp (Z2/E)1 2
For more details on aa--Gaussian distributions we refer to (Vershynin, 2018), Section 2.5.
We also need to introduce the space L LR Po) the space of essentially bounded functions, with finite || Ilxo norm, defined as:
|l inf {C [o(]f())] S CC)= =1
Note that the likelihood ratio l= dea for model (5) belongs to this space as soon as the tails of decrease faster than a standard Gaussian. For instance all the distributions with bounded support are included.
We often refer to the notion of events that happen with high probability, defined as follows.
Definition 7 (Events happening with high probability events). The sequence of events {Ea)den happens with high probability for every k 2 O, there exists de such that for d. 2 ds
1
P(Ea) < d
B.3. Hermite polynomials
In this section we group some facts about Hermite polynomials that will be needed for the subsequent proofs. For a comprehensive overview on the topic we refer to (McCullagh, 2018).
14

---

Feature learning from non-Gaussian inputs
Definition 8 (Hermite expansion). Consider function f IR + R that is square integrable with respect to the standard normal distribution p(x) = (1/V2T)) E e- Then there exists a unique sequence of real numbers {cx keN called Hermite coefficients such that
f(x) = Σ
Ck hx(2) K!
and cx(z) == : -=-=-woo.nl((())a(c))
where hú is the i-th probabilist's Hermite polynomial.
Definition 9 (Information exponent). Given any function f for which the Hermite expansion exists, its information exponent k* = (f) is the smallest index k 2 1 such that eg. +
Lemma 10 (Expectations of Triple and Quadruple products), For any i, k N, we have
i!j!k!
(i. R (ow..
(1 k
k
2
the denominator exists, that is if i k is even and the sum of any two of i,j and is strictly less than the third, or zero otherwise. What's 's more, fix the indices m1, m2, m2, m. N and order them such that my 2 m2 my 2 m44
Define M := (m1 + m2 + my + ma)/2. Then, given min(m3, ma) and (m11, ma, ma) = (e, Mamg (2))kmg (x) (1)) (2)], we have
m
(m3 + 2w))!m ! ma!! m4
T (m1, m2, ma, m4)=
Σ
(M mi m4 v))(M my MMM m2 v)M - vi)(mms vilu! ==
if .M 2N and the denominator exists. Otherwise, 12 (m1, m2: m3, m4) 0.
Lemma 11. Consider W1, W2 E g and set W1 U2- Then, for any i.j N
1x [h,(w1 z)hy (w> =)] a² 6ig.
To extend this orthogonality property to triple an quadruple products of Hermite polynomials, we use the same technique as 2015)
in Lemma 12 (Orthogonality property). Fix i, k N. Set 13 B, (4 (wi x 2:)e (w2 x 2))x (wz x)], Consider three different spikes wi, w2, wz E gd-1 and set p = wi W2, T == wi W3 and = w2 uzz Then, we have
ijjk!
I if the denominator exists, and 13- 0 otherwise.
p² (n
1= (j- 0)12 i
k
i
k - j + 28 j
2
2
2
Proof. We start by writing w2 sum of its projection along the first spike wi and its part in the orthogonal space, that is W2 = 1 + VI- puu, where u E Sd-1 and u = 0. Then, clearly, Wa x = pbı + v - p²ya. where := W1
and yu : u T.
We can do the same for the third spike, that is W/ = rw1 + /1 - 72u', where u' E S and u' ww = 0. Define you = x. We can now exploit the fact that any Hermite polynomial of sum can be written as the sum of Hermite polynomials the following way, that is for any E N we have that hm (azi+ /1 a²22) Em (""a" am)mmm (ee(z) (m for any [0. 1]. By applying this formula to Wa and W3, get
1/ = Σ
k.k
1
(2)))))=
of (VI 23) ( ( 73)
Ehi Mh (v1) (v1) ((w))-
E--==
(1)
(2)
where we have used the fact that is independent from Yv. and yor since wh and Lemma 11, we have
W2 n = O. Then, in view of Lemma 10
(1)
and
(2) = (j - -))((u
1
1
1 -1)
2
2
15

---

Feature learning from non-Gaussian inputs
Hence, imposing that L' =
obtain
13
IV P V1 VI-+² (u. (u-''))-=
+ 1
1:
20
8
t= ( 8)1
2
2
when the denominator exists and it vanishes otherwise. u u'
Therefore, the thesis follows from the fact that we have
TP
VI -1--211-P
Corollary 13. Assume that E is probability distribution and that G E (RR is even. Then, if w1, w2. wy E gd-1 and W1 # W2 # Wy, we obtain that
E. ("m x))G (w2 x))(w3 x)] = (w) w22 -2-1)
Proof. By considering the Hermite expansion of G' and é, it clearly follows that
G &'G.d. c
(G x)) (w2 x))(w3 + x]] =
E ((
[hs(wo x))5(W2 x: (w)))
(10)
iljlk!
k=
where, in view of Lemma 12, when = w1 W2 T == wi W0 and 1) == W2 Wa it holds that
i!jik!
[[{+)]
++ (17 (p))))
(11)
( k
i + k
1
+ 2t
t=
(j f)12
2
and then
E (0.1.)
(G x)) x)((wg . x²] = Σ
Σ
coeff p'r"~jtt(n - =p)--
kj t=
where
G 9 ck
coeff
i j (j 1)!2
1
28 8
2
2
2
We can compute way this sum scales with the dimension by recalling that since & probability distribution, we have of 1. However, since G is even, of of 0. The smallest indices such that the sum does not vanish are then I 0 and j K 1. Then, it turns that the only term which survives in the when 1 0 is the with t kg Therefore, up to checking that the other terms in the sum are dominated by the addendum corresponding to
A 0, i.j = kgg and / kī 1, we have
[G'(w x)G x)<(13 z)] 011
(0,
Corollary 14 (used formulas with same argument). Fix k: = 1 and i,j E N. Consider Ta given by 12 in the case of Wy = W2- Then,
ilj
I =
(5+1-1): p.
1 (
1
1 -
J +
2 2
2
when the denominator exists and I3 = otherwise. Moreover, for any i,j, k. E N, consider 14 = (ha
(w) with wi, W2, u3 E 8d-1
Then,
k (0) )
4 4= --E-I
E=0
where 13 is defined Lemma 10 and the indeces i, J, and + El are ordered.
16

---

Feature learning from non-Gaussian inputs
Proof. The first part is obtained by repeating the proof of Lemma 12 with w2 pmd wi 0, and using the fact that (0, Ma)lds-e(u x)] = Óz-t- For the second part, it is sufficient to repeat the usual procedure by decomposing both w2 and wg such that W2 pwi + vI - p and Wy = Twi + vI (11-7ap 7zu for u' E gd-1 and wi W2 = Wy = 0. Then, the the- sis follows from the second part of Lemma 12.
for E gd-1 and
C. Measures of non-Gaussianity for ICA
Whether the population loss has to be (globally) maximized minimized with respect to the weight vector E Sd-1 depends on the data distribution and the choice of the contrast function G. Therefore, we briefly comment on this, by referring (Hyvärinen & Oja, 2000) and (Hyvärinen et al., 2001) for more details. What does maximising non- non-Gaussianity mean? Equivalently: what is suitable measure for non-Gaussinn non-Gausiianity A standard approach is to choose the negentropy as a measure of the oon-Gausiannty of the random vector x' ~ P. i.e.
(x) = H (I - (x),,
where Hxx) := [Ey]lo((z)] and TGauss Pg- We assume that x has zero mean and that its covariance matrix is the identity. The negentropy, being Kullback-Leible divergence, is always non-negative. However, to compute J one needs access the density of the random vector which is usually unknown and not easy to approximate. Hence, the classical approach is to approximate negentropy indirectly. In Section 5 of (Hyvärinen et al., 2001), the authors show that the negentropy can be approximated by
[(w)= (E [G(w x)] E (G( )))
where even and does not grow too fast because of reasons of numerical stability. Typically, G is one of the contrast
functions from Equation (3). Note that this optimization problem slightly different from Equation (1) (which the setting
considered in practice, mentioned in (Hyvärinen & Oja, 2000)). The reason why is convenient to remove the square clear as soon as the Hermite expansion Equation (14) is performed, together with the fact E [G( = cg c-:
[(w)= (f(w) - E E[G(IGuuss)))
-(866 1/6 +
- G E.ck. a + - - g66)
2
h!
G ch-
2
2h
+ )(k22)
k!
hence the information exponent of this loss function is doubled with respect to L worsening the sample complexity for the search phase (geometrically, this can be seen as the square "flattening" the loss around the origin).
The removal the square comes with the side effect that it is not sufficient to just minimize the loss, since the sign becomes relevant: G needs to me maximised in case of platykurtic data since the cumulants of the data are smaller than the gaussian ones, E(G(w =x]] (T < 0, and conversely needs to be minimized in case of leptokurtic data distributions.
C.1. Details on the noisy ICA model model
Here we expand the population loss (2) in the specific case where the inputs are distributed according the spike cumulant model (5). We can expand the likelihood ratio É in series of Hermite polynomials in the sense of Definition 8, and write
& Eg Σ
(lu x)
= x)
Note that due to the symmetry of cej+1 = for any j E N. Moreover, due to the whitening matrix, also ci 0; so, apart from of = 1, the smallest non-zero coefficient is oy. the one corresponding first non-trivial cumulant, i.e. the kurtosis.
17

---

Feature learning from non-Gaussian inputs
Using the orthogonality identity for Hermite polynomials (Lemmaa 11 in Appendix B.3), we get
[(w) = Ep, [G(w- ((v x))
20
dd
Brolhs 20),(00 x)1
ilj!
1220
c 4!
= 4 % +
+ olaa))
where a = w : v the only order parameter of the system. For this calculation, we have assumed that of o 0. which holds for the contrast functions in Equation (3). This implies that the inference problem of detecting the spike has information exponent k* 4, in the sense of Definition 8 in Appendix B.3.
Note also that, for this specific data distribution, l can be computed explicitly, that is
E(v) Ey V B
Bexp
1 + (-2 2
B
+
(12)
+
3
2
For more details regarding this computation, we refer to Szekely al. (2024).
D. FastICA
In Figure 5 we show numerically that the large fluctuations seen in Figure 2 (middle) hold for various values of e. We run FastCCA on the noisy ICA model with size batch d3+e, for e = 0.5, 0.8. We plot ay v w in the first four iterations, plus the overlap at initialisation.
We prove now Proposition 1, which guarantees that the second order term in the iteration of FastICA does not only speed up the convergence of the algorithm, but helps in escaping the search phase too: indeed, when is not present, the recovered signal after the first giant step scales the random weight initialisation, i.e. aj (v w,)² = O(1/d). However, when the iteration of FastICA is complete, the relevant direction is learned in one single step, if the size-batch is infinite and the dimension of the inputs diverge.
Proof of Proposition (Infinite batch size). We start by recalling that, up to normalisation, the first iteration of FastICA in the population limit, i.e. for a batch of infinite size, reads
W1 Eix Gw x)] - B[G" x)] wo.
where the weight vector initialisation wi uniformly drawn in the unit sphere in Rd.. Here, we have considered the complete iteration, with both the first order and the second order terms. Since assume we are in high dimensional regime, we have that the overlap between the initialisation and the spike that we want to recover scales as ag wo v (l/vdd) We want to prove that, after the first step, the spike has been totally learned when d - +00, namely
af = (wi v)² = 1 - o(1).
To do so, we can expand in the series of Hermite polynomials the squared overlap
with क1 = W1
(13)
at sufficiently high order. First of all, take into account the first and the second order terms f := x G' (wo x) and = G (wp x.) w) and compute their expectations with respect to the data distribution.
We compute [ by noticing that f E Rd can be written as the sum of a spherical gradient and a derivative along the direction of the initialisation, in the following sense: if L(w, x) == G(w - x), we have that
f (1a wowm V n. L.
+ wowb V.. v..((w..)
VanaL(w, womang
woo T V..L(uu x))
w=
w=
1m0000
(+)
(x)2
18

---

Feature learning from non-Gaussian inputs
ned 32
n=d""1,e=20
05
0.0
FastiCA steps
Figure 5. Performance of FastICA on the spiked cumulant model for various sample complexity regimes. We run the standard CastICA updates on data drawn from noisy ICA model (5), At cach step of the algorithm, we draw new batch of size d with e = 0.2, 0.5, 0.8. We plot the overlap between the planted Non-Gaussian direction and the estimate produced by FastCCA. Theorem confirms this picture for a general data model.
This splitting has to be done since need unit vectors to be allowed to apply the orthogonality formulas in Section B. By using the likelihood ratio trick, we get
[[(+)] =E
2 Ow
G(w x))
G(w 2)) 1==ma
llv x)
Ow
G ne(w-uy m))
EowN((0,a))
Σ
Σ j!
hy (wo - x))
Sw
i!
10=00
0
E
(0.0)
xw)
Sw
ilj!
w==w
=,=00
0
CK ck
Σ
cl (k 1}1 0
k-1
V.Lo Σ
c (k 1)! 0
k-1
8w
k!
0=
w=
where Viwg == 0 Wy is the projection of the spike in the space orthogonal to Wo- Moreover, recall that k- is the information exponent of the population loss (2), in the sense of Definition 9. In the case of the spiked model defined in Section 2.2, already know that k 4. We keep the terms up to the sixth-order of expansion. Hence,
cf
G c 0 +
E [({+)]] = Ulwo
ag
+ +(ag))
3!
5!
On the other hand, also have to take into account the component of the gradient in the direction of the initialization. For this second part, we get
d
[[(+)2]
=w E
L L(Xwo, x)) A=1
wo E [G'(wo <wo x)]
dd
=wo E [G' (wo x)) 21 (wo x) x)]
E [ha(wo 2) - x0h5(v- (x)
=wo
dly!
By using Corollary we can expand up to six order, obtaining
c
& cl
E [
2
G
0
4
+
00)
+(ag))
3!
4!
5!
6!
19

---

Feature learning from non-Gaussian inputs
In conclusion, summing up the contributions simply have
9 0g
c wo
d 3!
&
[
V mio
ag
(
00
00
+
+
5!
4!
5!
6!
Therefore, it remains compute Blg].
To do so, is sufficient to expand G" in series of Hermite polynomials to obtain
Ge o
c
& 0B
wo (
4! 0 - 66 +
+ olag)
wp
0
h!
6!
Now that we have expanded both E[f] and Eig) up to sixth order, after long but simple calculation, we can conclude by writing
(cf
ag + o(ag)
=
312
This computation has shown that if the first step could possibly be infinitely large, then FastCCA would perfectly learn the spike up contributions which vanish when +
We observe now that happens when the iteration of FastICA consists only in the gradient term. The different update would be given by wi = E[ G' (wo =)]. In this case, it is not needed to expand up to sixth-order, since we can easily see that of scales as bad as the random initialisation. In particular,
= a = ai + +o(a)) = =2
. (g)))220 + o(a))) | (cZ o(1)
We going prove Theorem which makes explicit the sample complexity required learn the spike in single step with large, although finite, amount of data points. Recall that, given samples, the FastICA iteration reads
1 n n
W.
x Gww 1..1x)
x²)
Σ
G" (w,-1 x0)),,-
n2
G' z²) and g == 1 E. x" G" x" w0)) The key ideas of this strategy are borrowed from the proof of Theorem 1 and 2 in Dandi al. (2024), where they consider teacher-student setup and Gaussian inputs. First of all, we compute the expectations of we look at how much data we need to obtain concentration for both the quantities, namely what is the scaling for such that I and || 2 are approximated, up to sufficiently small error, by their expectations. Then, we get the scaling for of by simply computing the ratio between E[TT]] 12 and E[ W 11
Proof of Theorem 2 (Finite batch-size). Define f= 1 -
n
Wi and || 1| After that,
From now on, we consider n = ((ak* for 8 [0, In order to compute the expectation of Ii we can use the calculation performed in the proof of Proposition to obtain
1
E m = u. Elf - gl =
(B|x Glo'"G + x²)]
ElG"(un x²)))) wo))wo)
n
G ct
x- + x(axk-1))
0
- 1)!
where the last inequality is obtained by following the steps of Proposition without explicitly using that k* 4. Indeed,
G. ---
w- E[] E ap + (1- -11 a0
-(00--1) and
E[g] che ao + +(aa -11).
20

---

Feature learning from non-Gaussian inputs
We now have compute ||| w1
[1
B Elllg 2 E g]. Since
Bilg|
E
E
(x* z² Gwwu x²) G' (wo x²)
n
E
n
(wo (x" x²) G' (wo x²) G' (wo x²))
+
E n
1
n
||~0|| 1 G' - =²)²
n1
nnn
1)
[B
B | G"wa 2)²
n²
and it is possible to do the same calculation for f|ll and Ef g]. we get that
n(n 1)
u |?]
+ |B[]] [1]]
2 [R[g] [B([))
n² [B||=|| (wo x) BlG" =)] 2 [x (o))o z)G (wo =))):
(14)
We can compute the scaling for the last three terms by using Lemma 12. On the other hand, the scaling for each of the first three terms is obtained by expanding fl and E[g] up to sixth-order and repeating the same calculation of Proposition I while avoiding the evaluation h- = 4. We get that, for any n, we have
||]|w]|| ||²]
n2nn 1)
(c2) )2
2(1 2(1'-1) -1)
n
(k'
1)!2"
+ d
Σ
E Ck
+ O(da +
Σ
Ck ck
+ +0(a0' - of' 2200"'o''66 6f'66+00) 011)
k!
k!
KE2N+1
KEZN
(o e2
c2-)
C & k!
2(k"
Olda
kK* 1)12
au
n²
(k* 1)!22 4D
k
G" ck
+
+0(at)
200 6 c + +
kE2N
Hence, start to keep into account the different choices of 8 [0, 2). Recall that O(d"'-o) We need to know that is the addendum which dominates in the previous formula, depending on k" By looking at the scaling of each addendum, one 2(k
can see that there are two possible regimes, which correspond to the two cases ap ((K--)) $ d/n. In particular, when the signal coming from the contribution of the first line in (14) dominates. More precisely, get that
B|]l |1 θ (
8 E (0, 2]
B|[lw113 (cL cf-cct-2
(k* 1)12 d + ola0 2(K -1))
-(---1)
8-0
Hence, have the scaling of the expectations of Ti and || 1l. is left to see how much data samples are required to approximate well these expectations with high probability, i.e. for which sample complexity æ, and ||w]|| || concentrate. Thanks to Assumptions 1, and 3, the proof is analogous the one in Section B.4 from (Dandi et al., 2024), leading to:
If v Ef all + +gg Elg w] = O (lg(()) Vn
m B|m]ll
and
wy 2| [|| w1 |1 < 111.f112 - E ||]f[| + |1la|1
E||l|| +2ff EEf gll log(n)* d 11122
= (
d log(n)
log(d)6
log(n)²
nyn _ +
nyd
n
21

---

o(1), Feature learning from non-Gaussian inputs
Negative results: We start by analysing what happens when å (0,2). 2]. We have with high probability that
| || 2 Ellbar
|1
wn [| - w| ||9| ( - -0
d log(n))
+
log(d))
log(1e)
log(re)*
nyn
nvd
d
11/21///
which implies that
(
d
|
=
n2
Moreover, with high probability,
|T| S [[[=]]] + [B[m]] =1 O
+0 log(n)
31 = 0
+0 (log(m) vn
dik d(k- 13/2
v'n
d(k"
We have split the reasoning in cases.
1 If E (0, 1). with high probability #1 =0
Case
and then
11//2
d
d
1 88
ai
O
W
d
Ω
dl
means that the overlap at step one is bounded by above by a quantity which vanishes when the dimension of the inputs goes infinity, which is definitely a bad news for learning high dimensions with FastICA. Note however that this case bridges smoothly the cases 8 0 and the case 8 1, which correspond d pi"--, i.e. the regime in which online SGD able to learn.
If = [1, 2). the reason why FastCCA does not learn is not an adverse ration between two random, although concen- trated, quantities. In this case, the number of samples is not enough to reach concentration for πT, that is Therfore,
Case 2:
O [lg(m)
polylog(n)
( (polylo(d)
of
W1
The latter negative result is qualitatively different from the previous one, since it due the lack of concentration of Ti. The spike is not recovered at all, meaning that the scaling of 0ı bad as the one provided by the random weights initialisation.
Positive result: Let's consider 5 = 0, which means that we employee an extensive number of data points, namely ((dh)) ). Then,
(c.- )2 (k*
|m 2
[E[| W1 || + [B[[|]]] [13] [| 11
2(k--1)
a(a 1)
1)12 0g
On the other hand, we also have the following lower bound:
G
c- ot
m Til
(k* - 1)!
90
-1
o(a x-
In conclusion, we get
(c). c
21k- o(a 2(k -11)
(k* 1) Qn
a1
[1
(og c }2
2(k=
o(a 2(k
2(1- -1)- 1))
On
+
(k+
1)!2
meaning that this regime, when +00, the relevant direction perfectly recovered by the first iteration of CastICA.
22

---

Feature learning from non-Gaussian inputs
Overlap = d
0.35
smoothing
0.04
Loss ((a)
0.0175
0.25
0.0150
0.03
0.0125
0.15
0.02
0.0100
0.10
0.0075
0.01
0.0050
0.00
-0.5
0.0
0.5
1.0
-1.0
-0.5
0.0
0.5
-1.0
-0.5
Overlap
Overlap a
Overlap a
Figure 6. Benefit of smoothing the loss depending on different values of the parameter x: We show the smoothed loss Lx|G [G x]] for different scalings of the smoothing parameter A. We see the loss is most flat around the origin without smoothing, while increasing the smoothing parameter makes the loss more peaked around the origin, making it easier to escape the origin. See Appendix E for details.
E. Smoothing the landscape
For implementation purposes, we start by writing the formula for the spherical gradient of the smoothed loss, that
V (=)) = VanLi(G(w z)) = Esg [o' 10 Xa PP
(m x)) | V+ )
V
VT+
where is one of the components, for example the first one, of random vector ~ 2 This formula can be obtained by following the proof of Lemma 9 in (Damian et al., 2023). Note that we need to compute the spherical gradient in order to write the iteration of smoothed online SGD.
We now briefly comment on the benefit smoothing the loss landscape, which is given by taking into account the expectation of the contrast function evaluated in the dot product between the (normalized) vector w + Az and the vector of the inputs x, x. The smoothing operator is defined in Definition 3. Recall that 2 is vector of
instead of evaluating the loss in white noise and note that in Theorem 4 the parameter À is large, namely it scales with the dimension d of the inputs. The key intuition behind the smoothing operator is that, thanks to large X, it allows to evaluate the loss function in regions that are far from the updated weights w/. collecting non-local signal that alleviates the flatness of the saddle of the loss function at the beginning of learning, reducing the length of the search phase. In Figure 6 we plot the smoothed loss function, in the cases of X = ///4 A = dz/8 and 0. The first case corresponds to the optimal value of X, that is, the one for which the optimal thresholds of Theorem are obtained. The latter case corresponds to the absence of smoothing. can be interesting to recall that for À >> dz/4 we have the "noise" which dominates the "signal', in the sense of Equation (15), where the admissible values of X are chosen. For the plots in Figure we have used n = 10 number of samples, d = 10 and G(s) = -ha(s).
In the following section we will go through the rigorous argument that formalizes the heuristic derivation presented in Section 1.
E.1. Proof
We prove now Theorem 4. which provides the timescales for recovering the relevant direction with smoothed online SGD, Since the Gaussian version of this theorem has already been proved in (Damian et al., 2023), we will only emphasise the main differences between the Gaussian and non-Gane non-Gaussian case, which are mainly due to the presence of the (bounded) likelihood ratio and the fact that the analysis cannot take into account only the information exponent ki of the contrast function, but also the information exponent of the population loss, and those in general are not equal. However, we will conclude that smoothed online SGD obtains the optimal sample complexity only in the matched case, which is for k kg.
Proof of Theorem 4 (Escaping mediocrity). The proof is made of two parts. The first part is dedicated the computation of the scaling for the two main quantities that take part in the iteration of smoothed online SGD, i.e. the signal and the noise. The first part ends when we have collected all the ingredients to prove an analogue of the thesis of Lemma 15 in (Damian et al., 2023), which for us is given by (17). This formula guarantees that the overlap between the spike and the updated weight vector given by smoothed online SGD is precisely the sum the previous overlap, the signal, the noise and a suitable martingale term. The second part is dedicated to the application of standard arguments in probability theory (for the case of vanilla and smoothed SGD see e.g. (Ben Arous et al., 2021) and (Damian et al., 2023) respectively). from which
23

---

Feature learning from non-Gaussian inputs
is straight forward conclude that at some point, i.e. for a certain number of used data points, SGD reaches regime in which the signal dominates the noise and the spike is recovered.
Signal: We start by studying the scaling of the population (spherical) gradient V sph (w)). This scaling will depend on the information exponent kf of the population loss. To do so, consider the expansion as series of Hermite polynomials of the population loss
[(w))
C h!
ck
a
where kj is its information exponent. In Lemma 8 of (Damian al., 2023) it is shown that apply the smoothing operator to each term and obtain that, for any k 2 0,
a
a² 2 1² d:
Ex(a*) = sx(a, with sx(a, (1 + A2)b/2
1//2
a² and is even,
11/2
a² 2² and l: is odd.
By computation, we can obtain that the population (spherical) gradient is V Lx = (v -wwca (a), with
E d b! da Lx(aa)
(a) =
The scalar product of this population gradient term and the relevant direction, i.e. v. Ti (1 a²)ex(a), is called signal. It can be proved, like in Lemma 10 of (Damian et al., 2023), that
Skg - (a,
cx(a) RJ
v1 X X2
Noise: Now we take look at the noise term. Fix v gd-1, with u 1 w. Recall that we are using the notation 9 V (w, z). Since we have assumed that the likelihood ratio é is bounded (e.g. the case of inputs drawn according to the noisy ICA model (5)), we get
[uu g) Bx-A((.11)) [[u. g) x)] S |1e|loo R [(u g)
(15)
-(k2
min (1 X2, vd)
{1+
(16)
1 + X2
where we have used the fact that À [1, d/4] and Corollary in (Damian al., 2023). In the last inequality, thanks the choice À E [1, dd/4), we using the fact that
min (1+X², vd -1)-1) = + x2) (R] (1)
Note that we have obtain an explicit dependence on the information exponent kjj of the contrast function. We can estimate now the second moment of the smoothed spherical gradient, i.e. |[]gg]?], which we call noise. To do so, recall that g is a spherical gradient and then course g 1. w. To get an upper bound on the scaling of the noise, sufficient to fix an orthonormal basis in Id such that e = w. Then,
d
E||g| = Σ E [e: 9) S do + =E))
Now that we have the desired upper bound, it is possible to compute the scaling for the expectation of some p- moments of llgil,, for p-sufficiently large. This corresponds to Corollary 2 in (Damian et al., 2023).
By considering both the scaling of the signal and the noise, we can now state the Non-Gaussian version of Lemma 15 in W: + T9t we obtain that with high probability
(Damian et al., 2023), which tells us that if 04+1 == w Wht1 1 and we+1 = [w + n9ell
0++1 04 + (1 (ax) + 2, + (E( [1
(17)
24

---

Feature learning from non-Gaussian inputs
where (Ee)EN is martingale with mean defined
Ze == 17 (g: v) (gt U)] + 1r E8 [rd]
and
VI 17 19 -1)
Ft = [as + v (ge w)]
Hence, each step of smoothed online SGD, the updated overlap between the new weight vector and the spike is nothing but the previous overlap plus martingale term, the signal - whose scaling depends on kī and the noise whose scaling depends on k5- We have reached the setting of Lemma 16 of (Damian et al., 2023). From here on, the proof can be concluded by applying the argument that translates the signal and noise estimated into sample complexities, introduced in the setting of vanilla online SGD by (Ben Arous al., 2021). The proofs of Lemma 16 and Lemma 17 from (Damian al., 2023) can be repeated straight forwardly, leading the required sample complexities for recovering the spike.
F. A Additional numeries
n 2d
2d
n = 2d
n = 2d
-0.7
Training
loss
Test
-0.8
0
10
0
5
10
5
10
0
10
1.0
1/vd
+ ICA
overlap
0.5
0.0
10
0
5
10
5 steps
10
0
10
steps
steps
steps
Figure 7. Performance of FastICA noisy ICA model with Rademacher prior the latent without resampling the training data each step. Parameters: input dimension 25, signal-to-noise ratio B = 5, online learning.
25

---

Feature learning from non-Gaussian inputs
n = 2d
n 2d
n = 2d
n 2d
-0.700
-0.725-
loss
-0.550-
Train
-0.775-
Test
5
10
0
5
10
0
5
10
0
5
10
1.0
1/v d
overlap
ICA
0.5
PCA
5
10
0
5 steps
10
0
5 steps
10
0
5 steps
10
steps
Figure 8. Performance of online FastICA noisy ICA model with the "3not4" prior the latent variable. Under this prior, the latent has mean zero, unit variance, a non-trivial third-order cumulant, but zero fourth order cumulant; see Appendix A.5 for details. Parameters: input dimension 25, signal to-noise ratio online learning.
n= 2d
n= 2d
= 2d3
= 2d'
-0.70
H
loss
-0.22
Training
Test
-0.74
5
10
5
10
5
10
5
10
1.0
1/vd
overlap
+ ICA
0.5
0.0
5
10
5
10
5
10
0
5 steps
10
steps
steps
steps
Figure Performance online FastICA on noisy ICA model with Laplace prior on latent. Parameters: input dimension d = 25, signal-to-noise ratio B 5, online learning.
26


Reviews:
[{"review_id": "ZQZwjeXB7S", "reviewer": "Reviewer_oLx5", "summary": "Motivated by empirical observations that features learned by deep convolutional networks resemble those recovered by independent components analysis (ICA), this paper presents a concrete algorithmic sample complexity bound for various algorithms for recovering a non-Gaussian direction from $d$-dimensional data. Notably, the paper establishes a sample complexity bound of $n \\gtrsim d^4$ for a popular ICA algorithm (FastICA), and demonstrates that SGD in fact outperforms this algorithm on this task, attaining the information-theoretic lower bound of $n \\approx d^2$ (albeit through a smoothed loss).", "questions": "", "limitations": "", "rating": 4, "confidence": null, "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": "Thank you for your accurate comments and for the attention to the supplementary material, including the strategies of the proofs. Your suggestions offered valuable food for thought.\n\n>  I don't see any glaring weaknesses to the paper, beyond the possible restrictiveness of the single non-Gaussian direction model -- however, this is rather trite, as algorithmic results for other multi-index settings are yet generally ill-understood.\n\nThanks for your observation. We decided to focus on single-index models in this paper because of the prevalence of Gaussian single index studies in the recent past (Ben Arous et al. JMLR '21, Bietti et al. NeurIPS '22, Damian et al. NeurIPS '23, Damian et al. arXiv:2403.05529, Bardone & Goldt ICML '24, Wang et al. NeurIPS '24). We are considering on an extension to multi-spike models for our future work.\n\n> The authors should remember to put in the Impact Statement.\n\nThank you for the reminder, we will add it.\n\n> Something that might be helpful is a discussion of why identifying non-Gaussian directions might be desirable or might happen automatically. A concrete mathematical statement here might be to show recovering these directions improves generalization error.\n\nThe reviewer raises a really interesting point, which is why it is desirable to learn non-Gaussian directions in the first place. Fascinating work in theoretical neuroscience has established that natural images have a highly non-Gaussian statistical structure: while pixel intensities themselves may follow roughly Gaussian distributions, the relationships between pixels are strongly non-Gaussian. Specifically, natural images are sparse in certain bases: edges, contours, and textures are more prevalent than random pixel noise. Gabor-like filters, i.e. non-Gaussian directions, are efficient at capturing these features that yield sparse image representations, see for example reviews such as Simoncelli & Olshausen, Annu. Rev. Neurosci. 24:1193–216 (2001) or the book by Hyvärinen, Hurri and Hoyer (2009).\n\nMore recently, there have been a series of works investigating the importance of non-Gaussian input structures in machine learning from a mathematical perspective, showing that neural networks will learn them if they are relevant for the task, i.e. if they improve generalisation; see in particular Ingrosso & Goldt, PNAS '22; Bardone & Goldt, ICML '24; Lufkin et al. NeurIPS '24. \n\nWe will add an extended discussion of this issue to the revised version of the paper, should it be accepted.", "claims_and_evidence": "The theoretical claims in this paper are very well-supported. Various lower bounds and fundamental limits from prior work are outlined. The upper bounds in the paper are accompanied by simple numerical experiments demonstrating that they capture the correct scaling.", "methods_and_evaluation_criteria": "The numerical experiments are helpful for contextualizing the theory, and are documented in the appendix.", "theoretical_claims": "I did not check the proofs entirely. However, I checked the proof strategies of the main Theorems and they make sense to me. The numerical experiments also are helpful to sanity check the correctness of the theorems (such as demonstrating the necessity of smoothing the loss for SGD).", "experimental_designs_or_analyses": "A sufficient description of experiment set-ups is contained in the main paper, with additional details in Appendix A.", "relation_to_broader_scientific_literature": "This work belongs to the general category of literature concerned with understanding the statistical properties of machine learning problems and algorithms from the lens of \"feature learning\", i.e. how algorithms actually recover predictive features that, e.g., go beyond the linear regime (captured by the information exponent). This work in particular has potential value to the community in two main ways: 1. as stated in the paper, empirical observations have shown particular structure learned by CNNs that seem to match those given by earlier ICA algorithms--the results in this paper make this connection explicit, 2. the sample complexity bounds in this paper seem to properly close the statistical-to-computational gap for the non-Gaussian direction recovery problem, showing that both FastICA and vanilla online SGD are suboptimal, and that smoothed-loss SGD closes the gap.", "essential_references_not_discussed": "Not that I'm aware of.", "other_strengths_and_weaknesses": "I think this paper is written well and quite straightforward to digest, despite the technicality of the theoretical tools. The main insights and arguments are likely of interest to the feature learning theory community. I have described some of the strengths that stood out to me earlier. I don't see any glaring weaknesses to the paper, beyond the possible restrictiveness of the single non-Gaussian direction model--however, this is rather trite, as algorithmic results for other multi-index settings are yet generally ill-understood.", "other_comments_or_suggestions": "- The authors should remember to put in the Impact Statement.\n\n- This paper shows that ICA is computationally and statistically tractable. Something that might be helpful is a discussion of why identifying non-Gaussian directions might be desirable or might happen automatically, since it seems the ultimate goal is to demonstrate that neural networks might implicitly learn non-Gaussian directions. A concrete mathematical statement here might be to show recovering these directions improves generalization error, analogous to how weak recovery improves generalization in standard Gaussian single-index models. This is likely a hard problem, but food for thought.", "questions_for_authors": "None outstanding.", "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": "", "extra_scores": {}}, {"review_id": "tbTroUOdbD", "reviewer": "Reviewer_Jmtu", "summary": "The study quantifies the sample complexity of two learning algorithms: FastICA, Stochastic Gradient Descent (SGD). The key results are the following: FastICA requires at least $n \\gtrsim d^4$ samples to recover a single non-Gaussian direction in high-dimensional inputs. SGD outperforms FastICA in feature learning, achieving better sample complexity, particularly when smoothing techniques are used ($n \\gtrsim d^2$). In the real-world dataset ImageNet, the strong non-Gaussianity of image data helps mitigate FastICA’s inefficiency.\n\n## update after rebuttal\nThank you for the authors' response. The authors have addressed my main concerns. After reviewing the other reviewers' comments, I will keep my score and recommend accepting this paper.", "questions": "", "limitations": "", "rating": 4, "confidence": null, "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": "Thank you for your detailed feedback and your useful comments, which spurred us to run an additional experiment to link ICA with deep CNNs, and to provide two additional plots to provide intuition on the effect of smoothing, and on the intermediate regime of FastICA. We start with these points before addressing your remaining comments below. We hope our reply alleviates any remaining concern; if not, please let us know, otherwise we would really appreciate it if you increased your score. Thank you for your time!\n\n> Readability of the introduction & explanations on why Independent Component Analysis learns similar filters as deep CNNs in Fig. 1.\n\nThank you for this suggestion -- we now realise the introduction jumped from deep CNNs to ICA too quickly. To make the connection more clearly, we have conducted an **additional experiment**: we trained three deep CNNs (AlexNet, Resnet18 and DenseNet121) on ImageNet and computed the excess kurtosis of the dot products $s$ between first-layer convolutional filters and ImageNet patches. We found that the excess kurtosis of $s$, which corresponds to the objective function of ICA, sharply increases after about 1000 steps of SGD, precisely when Gabor filters form in the first layer. This demonstrates empirically that neural networks seek out non-Gaussian directions when Gabor filters form, akin to what happens when running ICA. We summarise our results in a revised Figure 1 in https://figshare.com/s/0e72e797306c1a3f216a, and that we hope will clarify the motivation for studying ICA.\n\n> The explanation of why SGD benefits from smoothing the contrast function could be clearer [...]\n\nThe key intuition behind the smoothing operator is that, thanks to large $\\lambda$ (than scales with the input dimension), it allows to evaluate the loss function in regions that are far from the iterate weight $w_t$, collecting non-local signal that alleviates the flatness of the saddle around $\\alpha=0$ of $\\mathcal L$, reducing the length of the search phase.  A new plot, Fig 3 of https://figshare.com/s/0e72e797306c1a3f216a, shows the reduced flatness of the smoothed loss at the origin.\n\n> Choice of $\\delta$ in the batch sizes $n=d^2, d^{3 + \\delta}$, and $d^4$\n\nWe could have chosen any $\\delta \\in (0,1)$, since Theorem 2 implies that the 'intermediate' regime shown in Fig. 2 (middle) holds for a batch size of $d^2 \\ll n \\ll d^3$. We will add Fig 2 of https://figshare.com/s/0e72e797306c1a3f216a which shows the behaviour for  $\\delta = 0.5, 0.8$.\n\n> The Hermite expansion, square integrability and contrast functions.\n\nWe merely need the contrast functions to be square-integrable with respect to the standard normal distribution $P_0$, so we need that the growth of $G$ is less than exponential. In general, we don't require they are integrable with respect to the Lebesgue measure on the real line (e.g., excess kurtosis is not Lebesgue integrable).\n\n> In Figure 3, the learning rate $\\eta$ is defined based on $k^*_1$, $k^*_2$, but how to tune it in practice?\n\nOur theorem 4 offers a way to make a reasonable guess for $\\eta$ as follows: $k^*_2$ is the information exponent of the contrast function and is therefore known. Meanwhile, $k^*_1$ depends on the likelihood and hence on the data distribution. However, we have the bound $k^*_1 \\geq k^*_2$. Furthermore, if we assume that data has a non-trivial fourth-order cumulant (a mild assumption in practice), we can set $k^*_1=4$ and hence our theorem gives a recommendation on how to scale the learning rate in practice.\n\n> The analysis is focused on learning one non-Gaussian feature. A discussion on how multiple features interact would be valuable.\n\nWe decided to focus on single-index models in this paper because of the prevalence of Gaussian single index studies in the recent past (Ben Arous et al. (2021), Damian et al. (2023), Damian et al. (2024), Bardone & Goldt (2024)). We are considering on an extension to multi-spike models for our future work.\n\n> Ensuring present tense consistency\n\nThanks, we fixed this.\n\n> Highlighting the generality of non-Gaussian distributions (Cramér's decomposition theorem) [...]\n\nYes, we agree that Cramér's theorem emphasizes the importance of studying non-Gaussianities, we will discuss it. \n\n> Additional empirical justification to explain that \"the growth of the excess kurtosis might compensate the poor sample complexity of FastICA in practice\"\n\nWe will move the plot from the appendix to Fig. 4 with ImageNet results. We agree that it does not fully explain the behaviour of FastICA on ImageNet, but the growth of excess kurtosis with input dimension hints at important finite-size effects, which are out of the purview of our asymptotic theory.\n\n> Explicit emphasis in the title that this paper analyses sample complexity [...].\n\nThank you for this suggestion; we will consider clarifying the title should the paper be accepted. We have already made changes to clarify our objectives and motivations in the introduction (see below).", "claims_and_evidence": "Yes, for the most part, the claims are well-supported by theoretical derivations and empirical experiments. However, more empirical justification to explain that \"the growth of the excess kurtosis might compensate the poor sample complexity of FastICA in practice\" would be useful, currently it is not fully explored.", "methods_and_evaluation_criteria": "Yes, the methods and evaluation criteria are appropriate. The ImageNet experiment is a reasonable real-world test, although the authors could have explored other high-dimensional datasets for broader generalizability.", "theoretical_claims": "The proofs appear mathematically rigorous. I have some related questions, see “Questions For Authors”.", "experimental_designs_or_analyses": "The experimental setups are well-documented, specifying parameter choices, dataset details, and batch sizes.\nDifferent regimes for FastICA (Figure 2) and SGD (Figure 3) are clearly explained. The spiked cumulant model used for evaluation is well-motivated. One potential issue is the color scales of the filters in Figure 1(b) could be adjusted to better match those in Figure 1(a) for better visual comparability.", "relation_to_broader_scientific_literature": "The findings in this paper could be relevant for designing more efficient feature extraction methods in deep learning.", "essential_references_not_discussed": "Not founded.", "other_strengths_and_weaknesses": "Strength:  \n1. The paper is rigorously written with strong theoretical contributions, which advances the understanding of feature learning dynamics in high-dimensional settings.  \n2. The link between ICA and CNN feature learning is well-motivated and provides a fresh perspective on the emergence of structured filters.  \n3. The theoretical results are logically structured, featuring clear mathematical formulations and providing essential technical background information.\n\nWeaknesses: \n1. An explicit emphasis in the title that this is a paper analysing sample complexity might make it clearer to the reader.  \n2. The writing of the introduction can be improved to increase readability.  For example, the first paragraph of the introduction gives a somewhat disjointed impression. After discussing Gabor filters, it shifts to a discussion of non-Gaussianity. SGD is also introduced quite abruptly. The transitions between these topics are non-smooth, making it difficult for the reader to follow the logical flow. A similar issue arises in the section describing the contributions of the paper, where several key terms and expressions appear without prior introduction. As a result, it is challenging for the reader to grasp the main focus of the work upon first reading the introduction. Smoother transitions and clearer structuring would enhance readability and coherence.  \n3. The explanation of why SGD benefits from smoothing the contrast function could be clearer, especially for readers less familiar with statistical-to-computational gaps. Additional visualisations would be better.  \n4. The analysis is focused on learning one non-Gaussian feature. A discussion on how multiple features interact would be valuable.", "other_comments_or_suggestions": "1. It would be beneficial to provide additional explanations on why Independent Component Analysis learns similar filters as deep convolutional neural networks in Fig. 1. Specifically, more details are needed to illustrate how this similarity can be observed directly in the figure.  \n2. Titles should begin with capitalized initial letters.  \n3. \"We finally **investigated** whether FastICA also exhibits an extended search phase on real data.\" -> \"We finally **investigate** whether FastICA also exhibits an extended search phase on real data.\" to maintain present tense consistency. \n4. I think perhaps highlighting the generality of non-Gaussian distributions (e.g., using Cramér's decomposition theorem) would emphasise the importance of this work more significantly.", "questions_for_authors": "1. The batch sizes are described as $n=d^2$, $d^3+\\delta$, $d^4$, but it is unclear how $\\delta=0.2$ was chosen. Is there any reason to choose $\\delta=0.2$? Does it generalize to different values?\n2. The Hermite expansion assumes $f(x)$ is square-integrable. However, many ICA contrast functions are not bounded (e.g., excess kurtosis). Could the authors clarify whether the expansion holds for all practical contrast functions?  \n3. In Figure 3, the learning rate $\\eta$ is defined based on $k_1^*$ and $k_2^*$, but how to tune it in practice?", "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": "", "extra_scores": {}}, {"review_id": "0Jt5gmtjDC", "reviewer": "Reviewer_eH2y", "summary": "This paper investigates the unsupervised learning method ICA as a simplified framework for feature learning in (deep) CNNs. In particular, the authors derive sample complexity thresholds for escaping the search phase of FastICA and SGD, considering a toy model of a dataset sampled from an isotropic distribution perturbed by a single non-Gaussian direction. Their results highlight the poor performance of ICA compared to SGD, with the latter being able, in principle, to achieve the computational threshold.", "questions": "", "limitations": "", "rating": 4, "confidence": null, "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": "Thank you for your careful feedback and for having checked both the numerical experiments and the proofs. Your suggestions will definitely help to clarify the paper.\n\n> One potential weakness of this work is the relative simplicity of the model. However, I do not consider this a major weakness, as it is compensated by the originality of the results in relation to the existing literature, providing a baseline for investigating more complex settings.\n\nIndeed, our data model is a simple single-index model. We decided to focus on single-index models in this paper because of the prevalence of **Gaussian** single index studies in the recent past (Ben Arous et al. JMLR '21, Bietti et al. NeurIPS '22, Damian et al. NeurIPS '23, Damian et al. arXiv:2403.05529, Bardone & Goldt ICML '24, Wang et al. NeurIPS '24). Here, we provide an extension to the non-Gaussian case, which is closer to realistic data. We are considering on an extension to multi-spike models for our future work.\n\n> Inconsistent notation in Fig. 1.B is used instead of (b).\n\nThanks for noticing, we fixed this.\n\n> I suggest writing explicitly what quantity is represented in the legends of Fig.4.\n\nThank you, we have added a title to the legends to stress that these are the widths of the patches extracted from ImageNet images.", "claims_and_evidence": "Every claim is supported by convincing evidence.", "methods_and_evaluation_criteria": "Although the proposed analysis focuses on a relatively simple model to derive quantitative results, it provides valuable insights into the performance of two popular algorithms used in practical applications.", "theoretical_claims": "I checked the correctness of the proofs of Proposition 1, Theorem 2 and Theorem 4. I do not have any issues to discuss.", "experimental_designs_or_analyses": "I have reviewed the experimental details related to the figures in the paper. The explanations are clear, and I have no issues to discuss.", "relation_to_broader_scientific_literature": "The paper extends the findings of Auddy & Yuan (2024) on the optimal algorithmic threshold for ICA, rigorously analyzing the performance of the most widely used algorithms in practice and comparing them to the known optimal result. Moreover, their proofs extend recent methods developed in the context of supervised learning with Gaussian inputs to the unsupervised non-Gaussian setting.", "essential_references_not_discussed": "I am not aware of any relevant related work that has been omitted.", "other_strengths_and_weaknesses": "The presentation is overall clear and easy to comprehend. The assumptions are correctly stated and well-motivated, and the results are supported by numerical simulations. One potential weakness of this work is the relative simplicity of the model, which perturbs the isotropic Gaussian distribution along a single direction and is distant from realistic image datasets. However, I do not consider this a major weakness, as it is compensated by the originality of the results in relation to the existing literature, providing a baseline for investigating more complex settings.", "other_comments_or_suggestions": "- Inconsistent notation in Fig. 1: $\\mathbf{B}$ is used instead of $\\mathbf{(b)}$\n- I suggest writing explicitly what quantity is represented in the legends of Fig. 4", "questions_for_authors": "I do not have additional questions", "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": "", "extra_scores": {}}]

Author rebuttal:

Meta-review: This work is motivated by the link between independent component analysis (ICA) and CNN feature learning, and proposes to advance the understanding of feature learning through the simple model of ICA, which seeks to recover the most non-Gaussian projection. It compares FastICA, a standard ICA algorithm, and stochastic gradient descent (SGD), which is commonly used in the training of deep neural networks, in terms of sample complexity. The derived performance bounds imply that, at the computational threshold of quadratic sample complexity, FastICA does not have concentrated update to allow it escape effectively from random initialization, whereas SGD, when properly smoothed in a data-dependent way, is able to recover the signal.

As remarked by the reviewers, the study is presented in a carefully structured manuscript, containing meaningful theoretical results under well-motivated assumptions, supported by clearly explained experiments. The simplicity of the theoretical setup is compensated by the originality of the results. The reviewers have maintained their original positive ratings after the rebuttal, during which the authors have provided convincing replies to the reviewers’ comments.
Decision: Accept (spotlight poster)