A Random Matrix Approach to Neural Networks

Introduction

In terms of practical applications, our findings shed light on the already incompletely understood extreme learning machines which have proved extremely efficient in handling machine learning problems involving large to huge datasets (Huang et al., 2012; Cambria et al., 2015) at a computationally affordable cost. But our objective is also to pave to path to the understanding of more involved neural network structures, featuring notably multiple layers and some steps of learning by means of backpropagation of the error.

These findings provide new insights into the roles played by the activation function $\sigma(\cdot)$ and the random distribution of the entries of $W$ in random feature maps as well as by the ridge-regression parameter $\gamma$ in the neural network performance. We notably exhibit and prove some peculiar behaviors, such as the impossibility for the network to carry out elementary Gaussian mixture classification tasks, when either the activation function or the random weights distribution are ill chosen.

Besides, for the practitioner, the theoretical formulas retrieved in this work allow for a fast offline tuning of the aforementioned hyperparameters of the neural network, notably when $T$ is not too large compared to $p$. The graphical results provided in the course of the article were particularly obtained within a $100$- to $500$-fold gain in computation time between theory and simulations.

The remainder of the article is structured as follows: in Section 2, we introduce the mathematical model of the system under investigation. Our main results are then described and discussed in Section 3, the proofs of which are deferred to Section 5. Section 4 discusses our main findings. The article closes on concluding remarks on envisioned extensions of the present work in Section 6. The appendix provides some intermediary lemmas of constant use throughout the proof section.

Reproducibility: Python 3 codes used to produce the results of Section 4 are available at https://github.com/Zhenyu-LIAO/RMT4ELM

Notations: The norm $\|\cdot\|$ is understood as the Euclidean norm for vectors and the operator norm for matrices, while the norm $\|\cdot\|_{F}$ is the Frobenius norm for matrices. All vectors in the article are understood as column vectors.

System Model

From a neural network viewpoint, the $n$ neurons of the network are the virtual units operating the mapping $W_{i\cdot}x\mapsto\sigma(W_{i\cdot}x)$ ($W_{i\cdot}$ being the $i$-th row of $W$), for $1\leq i\leq n$. The neural network then operates in two phases: a training phase where the regression matrix $\beta$ is learned based on a known input-output dataset pair $(X,Y)$ and a testing phase where, for $\beta$ now fixed, the network operates on a new input dataset $\hat{X}$ with corresponding unknown output $\hat{Y}$.

where we defined $\Sigma\equiv\sigma(WX)$. This follows from differentiating the mean square error along $\beta$ to obtain $0=\gamma\beta+\frac{1}{T}\sum_{i=1}^{T}\sigma(Wx_{i})(\beta^{\sf T}\sigma(Wx_{i})-y_{i})^{\sf T}$, so that $(\frac{1}{T}\Sigma\Sigma^{\sf T}+\gamma I_{n})\beta=\frac{1}{T}\Sigma Y^{\sf T}$ which, along with $(\frac{1}{T}\Sigma\Sigma^{\sf T}+\gamma I_{n})^{-1}\Sigma=\Sigma(\frac{1}{T}\Sigma^{\sf T}\Sigma+\gamma I_{T})^{-1}$, gives the result.

the resolvent of $\frac{1}{T}\Sigma^{\sf T}\Sigma$. The matrix $Q$ naturally appears as a key quantity in the performance analysis of the neural network. Notably, the mean-square error $E_{\rm train}$ on the training dataset $X$ is given by

Under the growth rate assumptions on $n,p,T$ taken below, it shall appear that the random variable $E_{\rm train}$ concentrates around its mean, letting then appear ${\rm E}[Q^{2}]$ as a central object in the asymptotic evaluation of $E_{\rm train}$.

where $\hat{\Sigma}=\sigma(W\hat{X})$ and $\beta$ is the same as used in (1) (and thus only depends on $(X,Y)$ and $\gamma$). One of the key questions in the analysis of such an elementary neural network lies in the determination of $\gamma$ which minimizes $E_{\rm test}$ (and is thus said to have good generalization performance). Notably, small $\gamma$ values are known to reduce $E_{\rm train}$ but to induce the popular overfitting issue which generally increases $E_{\rm test}$, while large $\gamma$ values engender both large values for $E_{\rm train}$ and $E_{\rm test}$.

From a mathematical standpoint though, the study of $E_{\rm test}$ brings forward some technical difficulties that do not allow for as a simple treatment through the present concentration of measure methodology as the study of $E_{\rm train}$. Nonetheless, the analysis of $E_{\rm train}$ allows at least for heuristic approaches to become available, which we shall exploit to propose an asymptotic deterministic approximation for $E_{\rm test}$.

From a technical standpoint, we shall make the following set of assumptions on the mapping $x\mapsto\sigma(Wx)$.

Under the notations of Assumption 1, we have in particular $W_{ij}\sim\mathcal{N}(0,1)$ if $\varphi(t)=t$ and $W_{ij}\sim\mathcal{U}(-1,1)$ (the uniform distribution on $$) if$\varphi(t)=-1+2\frac{1}{\sqrt{2\pi}}\int_{t}^{\infty}e^{-x^{2}}dx$($\varphi$is here a$\sqrt{2/\pi}$-Lipschitz map).

We further need the following regularity condition on the function $\sigma$.

The function $\sigma$ is Lipschitz continuous with parameter $\lambda_{\sigma}$.

This assumption holds for many of the activation functions traditionally considered in neural networks, such as sigmoid functions, the rectified linear unit $\sigma(t)=\max(t,0)$, or the absolute value operator.

When considering the interesting case of simultaneously large data and random features (or neurons), we shall then make the following growth rate assumptions.

while $\gamma,\lambda_{\sigma},\lambda_{\varphi}>0$ and $d$ are kept constant. In addition,

Main Results

As a standard preliminary step in the asymptotic random matrix analysis of the expectation ${\rm E}[Q]$ of the resolvent $Q=(\frac{1}{T}\Sigma^{\sf T}\Sigma+\gamma I_{T})^{-1}$, a convergence of quadratic forms based on the row vectors of $\Sigma$ is necessary (see e.g., (Marc̆enko and Pastur, 1967; Silverstein and Bai, 1995)). Such results are usually obtained by exploiting the independence (or linear dependence) in the vector entries. This not being the case here, as the entries of the vector $\sigma(X^{\sf T}w)$ are in general not independent, we resort to a concentration of measure approach, as advocated in (El Karoui, 2009). The following lemma, stated here in a non-asymptotic random matrix regime (that is, without necessarily resorting to Assumption 3), and thus of independent interest, provides this concentration result. For this lemma, we need first to define the following key matrix

of size $T\times T$, where $w\sim\mathcal{N}_{\varphi}(0,I_{p})$.

for $t_{0}\equiv|\sigma(0)|+\lambda_{\varphi}\lambda_{\sigma}\|X\|\sqrt{\frac{p}{T}}$ and $C,c>0$ independent of all other parameters. In particular, under the additional Assumption 3,

Note that this lemma partially extends concentration of measure results involving quadratic forms, see e.g., (Rudelson et al., 2013, Theorem 1.1), to non-linear vectors.

With this result in place, the standard resolvent approaches of random matrix theory apply, providing our main theoretical finding as follows.

Let Assumptions 1–3 hold and define $\bar{Q}$ as

where $\delta$ is implicitly defined as the unique positive solution to $\delta=\frac{1}{T}\operatorname{tr}\Phi\bar{Q}$. Then, for all $\varepsilon>0$, there exists $c>0$ such that

As a corollary of Theorem 1 along with a concentration argument on $\frac{1}{T}\operatorname{tr}Q$, we have the following result on the spectral measure of $\frac{1}{T}\Sigma^{\sf T}\Sigma$, which may be seen as a non-linear extension of (Silverstein and Bai, 1995) for which $\sigma(t)=t$.

Let Assumptions 1–3 hold and, for $\lambda_{1},\ldots,\lambda_{T}$ the eigenvalues of $\frac{1}{T}\Sigma^{\sf T}\Sigma$, define $\mu_{n}=\frac{1}{T}\sum_{i=1}^{T}{\bm{\delta}}_{\lambda_{i}}$. Then, for every bounded continuous function $f$, with probability one

However, as shall be shown in Section 3.3, and contrary to empirical covariance matrix models of the type $P^{\sf T}W^{\sf T}WP$, $\Phi$ explicitly depends on the distribution of $W_{ij}$ (that is, beyond its first two moments). Thus, the aforementioned linearization of $\frac{1}{T}\Sigma^{\sf T}\Sigma$, and subsequently the deterministic equivalent for $\mu_{n}$, are not universal with respect to the distribution of zero-mean unit variance $W_{ij}$. This is in striking contrast to the many linear random matrix models studied to date which often exhibit such universal behaviors. This property too will have deep consequences in the performance of neural networks as shall be shown through Figure 3 in Section 4 for an example where inappropriate choices for the law of $W$ lead to network failure to fulfill the regression task.

For convenience in the following, letting $\delta$ and $\Phi$ be defined as in Theorem 1, we shall denote

Theorem 1 provides the central step in the evaluation of $E_{\rm train}$, for which not only ${\rm E}[Q]$ but also ${\rm E}[Q^{2}]$ needs be estimated. This last ingredient is provided in the following proposition.

As an immediate consequence of Proposition 1, we have the following result on the training mean-square error of single-layer random neural networks.

Let Assumptions 1–3 hold and $\bar{Q}$, $\Psi$ be defined as in Theorem 1 and (3). Then, for all $\varepsilon>0$,

Since $\bar{Q}$ and $\Phi$ share the same orthogonal eigenvector basis, it appears that $E_{\rm train}$ depends on the alignment between the right singular vectors of $Y$ and the eigenvectors of $\Phi$, with weighting coefficients

where we denoted $\lambda_{i}=\lambda_{i}(\Psi)$, $1\leq i\leq T$, the eigenvalues of $\Psi$ (which depend on $\gamma$ through $\lambda_{i}(\Psi)=\frac{n}{T(1+\delta)}\lambda_{i}(\Phi)$). If $\liminf_{n}n/T>1$, it is easily seen that $\delta\to 0$ as $\gamma\to 0$, in which case $E_{\rm train}\to 0$ almost surely. However, in the more interesting case in practice where $\limsup_{n}n/T<1$, $\delta\to\infty$ as $\gamma\to 0$ and $E_{\rm train}$ consequently does not have a simple limit (see Section 4.3 for more discussion on this aspect).

Theorem 3 is also reminiscent of applied random matrix works on empirical covariance matrix models, such as (Bai and Silverstein, 2007; Kammoun et al., 2009), then further emphasizing the strong connection between the non-linear matrix $\sigma(WX)$ and its linear counterpart $W\Phi^{\frac{1}{2}}$.

2 Testing performance

As previously mentioned, harnessing the asymptotic testing performance $E_{\rm test}$ seems, to the best of the authors’ knowledge, out of current reach with the sole concentration of measure arguments used for the proof of the previous main results. Nonetheless, if not fully effective, these arguments allow for an intuitive derivation of a deterministic equivalent for $E_{\rm test}$, which is strongly supported by simulation results. We provide this result below under the form of a yet unproven claim, a heuristic derivation of which is provided at the end of Section 5.

where $w\sim\mathcal{N}_{\varphi}(0,I_{p})$. In particular, $\Phi=\Phi_{XX}$ and $\Psi=\Psi_{XX}$.

With these notations in place, we are in position to state our claimed result.

Let Assumptions 1–2 hold and $\hat{X},\hat{Y}$ satisfy the same conditions as $X,Y$ in Assumption 3. Then, for all $\varepsilon>0$,

While not immediate at first sight, one can confirm (using notably the relation $\Psi\bar{Q}+\gamma\bar{Q}=I_{T}$) that, for $(\hat{X},\hat{Y})=(X,Y)$, $\bar{E}_{\rm train}=\bar{E}_{\rm test}$, as expected.

In order to evaluate practically the results of Theorem 3 and Conjecture 1, it is a first step to be capable of estimating the values of $\Phi_{AB}$ for various $\sigma(\cdot)$ activation functions of practical interest. Such results, which call for completely different mathematical tools (mostly based on integration tricks), are provided in the subsequent section.

The evaluation of (4) can be obtained through various integration tricks for a wide family of mappings $\varphi(\cdot)$ and activation functions $\sigma(\cdot)$. The most popular activation functions in neural networks are sigmoid functions, such as $\sigma(t)={\rm erf}(t)\equiv\frac{2}{\sqrt{\pi}}\int_{0}^{t}e^{-u^{2}}du$, as well as the so-called rectified linear unit (ReLU) defined by $\sigma(t)=\max(t,0)$ which has been recently popularized as a result of its robust behavior in deep neural networks. In physical artificial neural networks implemented using light projections, $\sigma(t)=|t|$ is the preferred choice. Note that all aforementioned functions are Lipschitz continuous and therefore in accordance with Assumption 2.

In anticipation of these likely generalizations, we provide in Table 1 the values of $\Phi_{ab}$ for $w\sim\mathcal{N}(0,I_{p})$ (i.e., for $\varphi(t)=t$) and for a set of functions $\sigma(\cdot)$ not necessarily satisfying Assumption 2. Denoting $\Phi\equiv\Phi(\sigma(t))$, it is interesting to remark that, since $\arccos(x)=-\arcsin(x)+\frac{\pi}{2}$, $\Phi(\max(t,0))=\Phi(\frac{1}{2}t)+\Phi(\frac{1}{2}|t|)$. Also, $[\Phi(\cos(t))+\Phi(\sin(t))]_{a,b}=\exp(-\frac{1}{2}\|a-b\|^{2})$, a result reminiscent of (Rahimi and Recht, 2007).It is in particular not difficult to prove, based on our framework, that, as $n/T\to\infty$, a random neural network composed of $n/2$ neurons with activation function $\sigma(t)=\cos(t)$ and $n/2$ neurons with activation function $\sigma(t)=\sin(t)$ implements a Gaussian difference kernel. Finally, note that $\Phi({\rm erf}(\kappa t))\to\Phi({\rm sign}(t))$ as $\kappa\to\infty$, inducing that the extension by continuity of ${\rm erf}(\kappa t)$ to ${\rm sign}(t)$ propagates to their associated kernels.

where we defined $(a^{2})\equiv[a_{1}^{2},\ldots,a_{p}^{2}]^{\sf T}$.

It is already interesting to remark that, while classical random matrix models exhibit a well-known universality property — in the sense that their limiting spectral distribution is independent of the moments (higher than two) of the entries of the involved random matrix, here $W$ —, for $\sigma(\cdot)$ a polynomial of order two, $\Phi$ and thus $\mu_{n}$ strongly depend on ${\rm E}[W_{ij}^{k}]$ for $k=3,4$. We shall see in Section 4 that this remark has troubling consequences. We will notably infer (and confirm via simulations) that the studied neural network may provably fail to fulfill a specific task if the $W_{ij}$ are Bernoulli with zero mean and unit variance but succeed with possibly high performance if the $W_{ij}$ are standard Gaussian (which is explained by the disappearance or not of the term $(a^{\sf T}b)^{2}$ and $(a^{2})^{\sf T}(b^{2})$ in (3.3) if $m_{4}=m_{2}^{2}$).

Practical Outcomes

We discuss in this section the outcomes of our main results in terms of neural network application. The technical discussions on Theorem 1 and Proposition 1 will be made in the course of their respective proofs in Section 5.

We first provide in this section a simulation corroborating the findings of Theorem 3 and suggesting the validity of Conjecture 1. To this end, we consider the task of classifying the popular MNIST image database (LeCun, Cortes and Burges, 1998), composed of grayscale handwritten digits of size $28\times 28$, with a neural network composed of $n=512$ units and standard Gaussian $W$. We represent here each image as a $p=784$-size vector; $1\,024$ images of sevens and $1\,024$ images of nines were extracted from the database and were evenly split in $512$ training and test images, respectively. The database images were jointly centered and scaled so to fall close to the setting of Assumption 3 on $X$ and $\hat{X}$ (an admissible preprocessing intervention). The columns of the output values $Y$ and $\hat{Y}$ were taken as unidimensional ($d=1$) with $Y_{1j},\hat{Y}_{1j}\in\{-1,1\}$ depending on the image class. Figure 1 displays the simulated (averaged over $100$ realizations of $W$) versus theoretical values of $E_{\rm train}$ and $E_{\rm test}$ for three choices of Lipschitz continuous functions $\sigma(\cdot)$, as a function of $\gamma$.

Note that a perfect match between theory and practice is observed, for both $E_{\rm train}$ and $E_{\rm test}$, which is a strong indicator of both the validity of Conjecture 1 and the adequacy of Assumption 3 to the MNIST dataset.

We subsequently provide in Figure 2 the comparison between theoretical formulas and practical simulations for a set of functions $\sigma(\cdot)$ which do not satisfy Assumption 2, i.e., either discontinuous or non-Lipschitz maps. The closeness between both sets of curves is again remarkably good, although to a lesser extent than for the Lipschitz continuous functions of Figure 1. Also, the achieved performances are generally worse than those observed in Figure 1.

It should be noted that the performance estimates provided by Theorem 3 and Conjecture 1 can be efficiently implemented at low computational cost in practice. Indeed, by diagonalizing $\Phi$ (which is a marginal cost independent of $\gamma$), $\bar{E}_{\rm train}$ can be computed for all $\gamma$ through mere vector operations; similarly $\bar{E}_{\rm test}$ is obtained by the marginal cost of a basis change of $\Phi_{\hat{X}X}$ and the matrix product $\Phi_{X\hat{X}}\Phi_{\hat{X}X}$, all remaining operations being accessible through vector operations. As a consequence, the simulation durations to generate the aforementioned theoretical curves using the linked Python script were found to be $100$ to $500$ times faster than to generate the simulated network performances. Beyond their theoretical interest, the provided formulas therefore allow for an efficient offline tuning of the network hyperparameters, notably the choice of an appropriate value for the ridge-regression parameter $\gamma$.

2 The underlying kernel

Theorem 1 and the subsequent theoretical findings importantly reveal that the neural network performances are directly related to the Gram matrix $\Phi$, which acts as a deterministic kernel on the dataset $X$. This is in fact a well-known result found e.g., in (Williams, 1998) where it is shown that, as $n\to\infty$ alone, the neural network behaves as a mere kernel operator (this observation is retrieved here in the subsequent Section 4.3). This remark was then put at an advantage in (Rahimi and Recht, 2007) and subsequent works, where random feature maps of the type $x\mapsto\sigma(Wx)$ are proposed as a computationally efficient proxy to evaluate kernels $(x,y)\mapsto\Phi(x,y)$.

As discussed previously, the formulas for $\bar{E}_{\rm train}$ and $\bar{E}_{\rm test}$ suggest that good performances are achieved if the dominant eigenvectors of $\Phi$ show a good alignment to $Y$ (and similarly for $\Phi_{X\hat{X}}$ and $\hat{Y}$). This naturally drives us to finding a priori simple regression tasks where ill-choices of $\Phi$ may annihilate the neural network performance. Following recent works on the asymptotic performance analysis of kernel methods for Gaussian mixture models (Couillet and Benaych-Georges, 2016; Zhenyu Liao, 2017; Mai and Couillet, 2017) and (Couillet and Kammoun, 2016), we describe here such a task.

Let $x_{1},\ldots,x_{T/2}\sim\mathcal{N}(0,\frac{1}{p}C_{1})$ and $x_{T/2+1},\ldots,x_{T}\sim\mathcal{N}(0,\frac{1}{p}C_{2})$ where $C_{1}$ and $C_{2}$ are such that $\operatorname{tr}C_{1}=\operatorname{tr}C_{2}$, $\|C_{1}\|,\|C_{2}\|$ are bounded, and $\operatorname{tr}(C_{1}-C_{2})^{2}=O(p)$. Accordingly, $y_{1},\ldots,y_{T/2+1}=-1$ and $y_{T/2+1},\ldots,y_{T}=1$. It is proved in the aforementioned articles that, under these conditions, it is theoretically possible, in the large $p,T$ limit, to classify the data using a kernel least-square support vector machine (that is, with a training dataset) or with a kernel spectral clustering method (that is, in a completely unsupervised manner) with a non-trivial limiting error probability (i.e., neither zero nor one). This scenario has the interesting feature that $x_{i}^{\sf T}x_{j}\to 0$ almost surely for all $i\neq j$ while $\|x_{i}\|^{2}-\frac{1}{p}\operatorname{tr}(\frac{1}{2}C_{1}+\frac{1}{2}C_{2})\to 0$, almost surely, irrespective of the class of $x_{i}$, thereby allowing for a Taylor expansion of the non-linear kernels as early proposed in (El Karoui, 2010).

where only the last three functions (only found in the expression of $\Phi_{ab}$ corresponding to $\sigma(t)=\max(t,0)$, $|t|$, or $\cos(t)$) exhibit a quadratic term.

More surprisingly maybe, recalling now Equation (3.3) which considers non-necessarily Gaussian $W_{ij}$ with moments $m_{k}$ of order $k$, a more refined analysis shows that the aforementioned Gaussian mixture classification task will fail if $m_{3}=0$ and $m_{4}=m_{2}^{2}$, so for instance for $W_{ij}\in\{-1,1\}$ Bernoulli with parameter $\frac{1}{2}$. The performance comparison of this scenario is shown in the top part of Figure 3 for $\sigma(t)=-\frac{1}{2}t^{2}+1$ and $C_{1}=\operatorname{\rm diag}(I_{p/2},4I_{p/2})$, $C_{2}=\operatorname{\rm diag}(4I_{p/2},I_{p/2})$, for $W_{ij}\sim\mathcal{N}(0,1)$ and $W_{ij}\sim{\rm Bern}$ (that is, Bernoulli $\{(-1,\frac{1}{2}),(1,\frac{1}{2})\}$). The choice of $\sigma(t)=\zeta_{2}t^{2}+\zeta_{1}t+\zeta_{0}$ with $\zeta_{1}=0$ is motivated by (Couillet and Benaych-Georges, 2016; Couillet and Kammoun, 2016) where it is shown, in a somewhat different setting, that this choice is optimal for class recovery. Note that, while the test performances are overall rather weak in this setting, for $W_{ij}\sim\mathcal{N}(0,1)$, $E_{\rm test}$ drops below one (the amplitude of the $\hat{Y}_{ij}$), thereby indicating that non-trivial classification is performed. This is not so for the Bernoulli $W_{ij}\sim{\rm Bern}$ case where $E_{\rm test}$ is systematically greater than |\hat{Y}_{ij}|{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}=1}. This is theoretically explained by the fact that, from Equation (3.3), $\Phi_{ij}$ contains structural information about the data classes through the term $2m_{2}^{2}(x_{i}^{\sf T}x_{j})^{2}+(m_{4}-3m_{2}^{2})(x_{i}^{2})^{\sf T}(x_{j}^{2})$ which induces an information-plus-noise model for $\Phi$ as long as $2m_{2}^{2}+(m_{4}-3m_{2}^{2})\neq 0$, i.e., $m_{4}\neq m_{2}^{2}$ (see (Couillet and Benaych-Georges, 2016) for details). This is visually seen in the bottom part of Figure 3 where the Gaussian scenario presents an isolated eigenvalue for $\Phi$ with corresponding structured eigenvector, which is not the case of the Bernoulli scenario. To complete this discussion, it appears relevant in the present setting to choose $W_{ij}$ in such a way that $m_{4}-m_{2}^{2}$ is far from zero, thus suggesting the interest of heavy-tailed distributions. To confirm this prediction, Figure 3 additionally displays the performance achieved and the spectrum of $\Phi$ observed for $W_{ij}\sim{\rm Stud}$, that is, following a Student-t distribution with degree of freedom $\nu=7$ normalized to unit variance (in this case $m_{2}=1$ and $m_{4}=5$). Figure 3 confirms the large superiority of this choice over the Gaussian case (note nonetheless the slight inaccuracy of our theoretical formulas in this case, which is likely due to too small values of $p,n,T$ to accommodate $W_{ij}$ with higher order moments, an observation which is confirmed in simulations when letting $\nu$ be even smaller).

3 Limiting cases

We have suggested that $\Phi$ contains, in its dominant eigenmodes, all the usable information describing $X$. In the Gaussian mixture example above, it was notably shown that $\Phi$ may completely fail to contain this information, resulting in the impossibility to perform a classification task, even if one were to take infinitely many neurons in the network. For $\Phi$ containing useful information about $X$, it is intuitive to expect that both $\inf_{\gamma}\bar{E}_{\rm train}$ and $\inf_{\gamma}\bar{E}_{\rm test}$ become smaller as $n/T$ and $n/p$ become large. It is in fact easy to see that, if $\Phi$ is invertible (which is likely to occur in most cases if $\liminf_{n}T/p>1$), then

and we fall back on the performance of a classical kernel regression. It is interesting in particular to note that, as the number of neurons $n$ becomes large, the effect of $\gamma$ on $E_{\rm test}$ flattens out. Therefore, a smart choice of $\gamma$ is only relevant for small (and thus computationally more efficient) neuron layers. This observation is depicted in Figure 4 where it is made clear that a growth of $n$ reduces $E_{\rm train}$ to zero while $E_{\rm test}$ saturates to a non-zero limit which becomes increasingly irrespective of $\gamma$. Note additionally the interesting phenomenon occurring for $n\leq T$ where too small values of $\gamma$ induce important performance losses, thereby suggesting a strong importance of proper choices of $\gamma$ in this regime.

A phase transition therefore exists whereby $\delta$ assumes a finite positive value in the small $\gamma$ limit if $r/n<1$, or scales like $1/\gamma$ otherwise.

If instead $r>n$ (which is the most likely outcome in practice), as $\gamma\to 0$, $\bar{Q}\sim\frac{1}{\gamma}(\frac{n}{T}\frac{\Phi}{\Delta}+I_{T})^{-1}$ and thus

where $\Psi_{\Delta}=\frac{n}{T}\frac{\Phi}{\Delta}$ and $Q_{\Delta}=(\frac{n}{T}\frac{\Phi}{\Delta}+I_{T})^{-1}$.

These results suggest that neural networks should be designed both in a way that reduces the rank of $\Phi$ while maintaining a strong alignment between the dominant eigenvectors of $\Phi$ and the output matrix $Y$.

Proof of the Main Results

In the remainder, we shall use extensively the following notations:

Finally, because of exchangeability, it shall often be convenient to work with the generic random vector $w\sim\mathcal{N}_{\varphi}(0,I_{T})$, the random vector $\sigma$ distributed as any of the $\sigma_{i}$’s, the random matrix $\Sigma_{-}$ distributed as any of the $\Sigma_{-i}$’s, and with the random matrix $Q_{-}$ distributed as any of the $Q_{-i}$’s.

where $C,c>0$ are independent of $d$ and $\lambda_{f}$. As a corollary (see e.g., (Ledoux, 2005, Proposition 1.10)), for every $k\geq 1$,

Notations: In all subsequent lemmas and proofs, the letters $c,c_{i},C,C_{i}>0$ will be used interchangeably as positive constants independent of the key equation parameters (notably $n$ and $t$ below) and may be reused from line to line. Additionally, the variable $\varepsilon>0$ will denote any small positive number; the variables $c,c_{i},C,C_{i}$ may depend on $\varepsilon$.

We start by recalling the first part of the statement of Lemma 1 and subsequently providing its proof.

for $t_{0}\equiv|\sigma(0)|+\lambda_{\varphi}\lambda_{\sigma}\|X\|\sqrt{\frac{p}{T}}$ and $C,c>0$ independent of all other parameters.

The layout of the proof is as follows: since the application $w\mapsto\frac{1}{T}\sigma^{\sf T}A\sigma$ is “quadratic” in $w$ and thus not Lipschitz (therefore not allowing for a natural transfer of the concentration of $w$ to $\frac{1}{T}\sigma^{\sf T}A\sigma$), we first prove that $\frac{1}{\sqrt{T}}\|\sigma\|$ satisfies a concentration inequality, which provides a high probability $O(1)$ bound on $\frac{1}{\sqrt{T}}\|\sigma\|$. Conditioning on this event, the map $w\mapsto\frac{1}{\sqrt{T}}\sigma^{\sf T}A\sigma$ can then be shown to be Lipschitz (by isolating one of the $\sigma$ terms for bounding and the other one for retrieving the Lipschitz character) and, up to an appropriate control of concentration results under conditioning, the result is obtained.

for some $c,C>0$ independent of all parameters.

Finally, using again the Lipschitz character of $\sigma(w^{\sf T}X)$,

which, with the remark $t\geq 4t_{0}\Rightarrow(t-t_{0})^{2}\geq t^{2}/2$, may be equivalently stated as

As a side (but important) remark, note that, since

and thus, since $\|\cdot\|_{F}\geq\|\cdot\|$, we have

Thus, in particular, under the additional Assumption 3, with high probability, the operator norm of $\frac{\Sigma}{\sqrt{T}}$ cannot exceed a rate $\sqrt{T}$.

The aforementioned control of $\|\Sigma\|$ arises from the bound $\|\Sigma\|\leq\|\Sigma\|_{F}$ which may be quite loose (by as much as a factor $\sqrt{T}$). Intuitively, under the supplementary Assumption 3, if ${\rm E}[\sigma]\neq 0$, then $\frac{\Sigma}{\sqrt{T}}$ is “dominated” by the matrix $\frac{1}{\sqrt{T}}{\rm E}[\sigma]1_{T}^{\sf T}$, the operator norm of which is indeed of order $\sqrt{n}$ and the bound is tight. If $\sigma(t)=t$ and ${\rm E}[W_{ij}]=0$, we however know that $\|\frac{\Sigma}{\sqrt{T}}\|=O(1)$ (Bai and Silverstein, 1998). One is tempted to believe that, more generally, if ${\rm E}[\sigma]=0$, then $\|\frac{\Sigma}{\sqrt{T}}\|$ should remain of this order. And, if instead ${\rm E}[\sigma]\neq 0$, the contribution of $\frac{1}{\sqrt{T}}{\rm E}[\sigma]1_{T}^{\sf T}$ should merely engender a single large amplitude isolate singular value in the spectrum of $\frac{\Sigma}{\sqrt{T}}$ and the other singular values remain of order $O(1)$. These intuitions are not captured by our concentration of measure approach.

Since $\Sigma=\sigma(WX)$ is an entry-wise operation, concentration results with respect to the Frobenius norm are natural, where with respect to the operator norm are hardly accessible.

Back to our present considerations, let us define the probability space $\mathcal{A}_{K}=\{w,~{}\|\sigma(w^{\sf T}X)\|\leq K\sqrt{T}\}$. Conditioning the random variable of interest in Lemma 2 with respect to $\mathcal{A}_{K}$ and its complementary $\mathcal{A}_{K}^{c}$, for some $K\geq 4t_{0}$, gives

so that, with the same remark as before, for $t\geq\frac{4\Delta}{KT}$,

To avoid the condition $t\geq\frac{4\Delta}{KT}$, we use the fact that, probabilities being lower than one, it suffices to replace $C$ by $\lambda C$ with $\lambda\geq 1$ such that

The above inequality holds if we take for instance $\lambda=\frac{1}{C}e^{\frac{18C^{2}}{c}}$ since then $t\leq\frac{4\Delta}{KT}\leq\frac{24C\lambda_{\varphi}^{2}\lambda_{\sigma}^{2}\|X\|^{2}}{cKT}\leq\frac{6C\lambda_{\varphi}\lambda_{\sigma}\|X\|}{c\sqrt{pT}}$ (using successively $\Delta\geq\frac{6C}{c}\lambda_{\varphi}^{2}\lambda_{\sigma}^{2}\|X\|^{2}$ and $K\geq 4\lambda_{\sigma}\lambda_{\varphi}\|X\|\sqrt{\frac{p}{T}}$) and thus

Therefore, setting $\lambda=\max(1,\frac{1}{C}e^{\frac{{C^{\prime}}^{2}c}{2}})$, we get for every $t>0$

which, together with the inequality $P(\mathcal{A}_{K}^{c})\leq Ce^{-\frac{cTK^{2}}{2\lambda_{\varphi}^{2}\lambda_{\sigma}^{2}\|X\|^{2}}}$, gives

Indeed, if $4t_{0}\leq\sqrt{t}$ then $\min(t^{2}/K^{2},K^{2})=t$, while if $4t_{0}\geq\sqrt{t}$ then $\min(t^{2}/K^{2},K^{2})=\min(t^{2}/16t_{0}^{2},16t_{0}^{2})=t^{2}/16t_{0}^{2}$. ∎

As a corollary of Lemma 2, we have the following control of the moments of $\frac{1}{T}\sigma^{\sf T}A\sigma$.

with $t_{0}=|\sigma(0)|+\lambda_{\sigma}\lambda_{\varphi}\|X\|\sqrt{\frac{p}{T}}$, $\eta=\|X\|\lambda_{\sigma}\lambda_{\varphi}$, and $C_{1},C_{2}>0$ independent of the other parameters. In particular, under the additional Assumption 3,

We use the fact that, for a nonnegative random variable $Y$, ${\rm E}[Y]=\int_{0}^{\infty}P(Y>t)dt$, so that

which, along with the boundedness of the integrals, concludes the proof. ∎

Beyond concentration results on functions of the vector $\sigma$, we also have the following convenient property for functions of the matrix $\Sigma$.

for some $C,c>0$. In particular, under the additional Assumption 3,

for some $C,c>0$. Let’s consider in particular $g:W\mapsto f(\Sigma/\sqrt{T})$ and remark that

We can apply Lemma 3 for $f:R\mapsto\frac{1}{T}\operatorname{tr}(R^{\sf T}R-zI_{T})^{-1}$, since we have

Lemma 3 also allows for an important application of Lemma 2 as follows.

for some $C,c>0$ independent of the other parameters.

Let $f:R\mapsto\frac{1}{T}\sigma^{\sf T}A(R^{\sf T}R+\gamma I_{T})^{-1}B\sigma$. Reproducing the proof of Corollary 2, conditionally to $\frac{1}{T}\|\sigma\|^{2}\leq K$ for any arbitrary large enough $K>0$, it appears that $f$ is Lipschitz with parameter of order $O(1)$. Along with (7) and Assumption 3, this thus ensures that

for some $C,c>0$. We may then apply Lemma 1 on the bounded norm matrix $A{\rm E}[Q_{-}]B$ to further find that

As a further corollary of Lemma 3, we have the following concentration result on the training mean-square error of the neural network under study.

for some $C,c>0$ independent of the other parameters.

We apply Lemma 3 to the mapping $f:R\mapsto\frac{1}{T}\operatorname{tr}Y^{\sf T}Y(R^{\sf T}R+\gamma I_{T})^{-2}$. Denoting $Q=(R^{\sf T}R+\gamma I_{T})^{-1}$ and $Q^{H}=((R+H)^{\sf T}(R+H)+\gamma I_{T})^{-1}$, remark indeed that

As $\|Q^{H}(R+H)^{\sf T}\|=\sqrt{\|Q^{H}(R+H)^{\sf T}(R+H)Q^{H}\|}$ and $\|RQ\|=\sqrt{\|QR^{\sf T}RQ\|}$ are bounded and $\frac{1}{T}\operatorname{tr}Y^{\sf T}Y$ is also bounded by Assumption 3, this implies

for some $C>0$. The function $f$ is thus Lipschitz with parameter independent of $n$, which allows us to conclude using Lemma 3. ∎

The aforementioned concentration results are the building blocks of the proofs of Theorem 1–3 which, under all Assumptions 1–3, are established using standard random matrix approaches.

2 Asymptotic Equivalents

This section is dedicated to a first characterization of ${\rm E}[Q]$, in the “simultaneously large” $n,p,T$ regime. This preliminary step is classical in studying resolvents in random matrix theory as the direct comparison of ${\rm E}[Q]$ to $\bar{Q}$ with the implicit $\delta$ may be cumbersome. To this end, let us thus define the intermediary deterministic matrix

with $\alpha\equiv\frac{1}{T}\operatorname{tr}\Phi{\rm E}[Q_{-}]$, where we recall that $Q_{-}$ is a random matrix distributed as, say, $(\frac{1}{T}\Sigma^{\sf T}\Sigma-\frac{1}{T}\sigma_{1}\sigma_{1}^{\sf T}+\gamma I_{T})^{-1}$.

First note that, since $\frac{1}{T}\operatorname{tr}\Phi={\rm E}[\frac{1}{T}\|\sigma\|^{2}]$ and, from (7) and Assumption 3, $P(\frac{1}{T}\|\sigma\|^{2}>t)\leq Ce^{-cnt^{2}}$ for all large $t$, we find that $\frac{1}{T}\operatorname{tr}\Phi=\int_{0}^{\infty}t^{2}P(\frac{1}{T}\|\sigma\|^{2}>t)dt\leq C^{\prime}$ for some constant $C^{\prime}$. Thus, $\alpha\leq\|{\rm E}[Q_{-}]\|\frac{1}{T}\operatorname{tr}\Phi\leq\frac{C^{\prime}}{\gamma}$ is uniformly bounded.

which, from Lemma 6, gives, for $Q_{-i}=(\frac{1}{T}\Sigma^{\sf T}\Sigma-\frac{1}{T}\sigma_{i}\sigma_{i}^{\sf T}+\gamma I_{T})^{-1}$,

Note now, from the independence of $Q_{-i}$ and $\sigma_{i}\sigma_{i}^{\sf T}$, that the second right-hand side expectation is simply ${\rm E}[Q_{-i}]\Phi$. Also, exploiting Lemma 6 in reverse on the rightmost term, this gives

We study the two right-hand side terms of (5.2.1) independently.

For the first term, since $Q-Q_{-i}=-Q\frac{1}{T}\sigma_{i}\sigma_{i}^{\sf T}Q_{-i}$,

where we used again Lemma 6 in reverse. Denoting $D=\operatorname{\rm diag}(\{1+\frac{1}{T}\sigma_{i}^{\sf T}Q_{-i}\sigma_{i}\}_{i=1}^{n})$, this can be compactly written

for some $C,c>0$, so in particular, recalling that $\alpha\leq C^{\prime}$ for some constant $C^{\prime}>0$,

As a consequence of all the above (and of the boundedness of $\alpha$), we have that, for some $c>0$,

where $D_{2}={\rm diag}(\{\frac{1}{T}\sigma_{i}^{\sf T}Q_{-i}\sigma_{i}-\alpha\}_{i=1}^{n})$. Of course, since we also have $-aa^{\sf T}-bb^{\sf T}\preceq ab^{\sf T}+ba^{\sf T}$ (from $(a+b)(a+b)^{\sf T}\succeq 0$), we have symmetrically

so that, with a similar reasoning as in the proof of Corollary 1,

where we additionally used $\|Q\Sigma\|\leq\sqrt{T}$ in the first inequality.

Together with (5.2.1), we thus conclude that

where the first equality holds by exchangeability arguments.

where |\frac{1}{T}\operatorname{tr}\Phi({\rm E}[Q_{-}]-{\rm E}[Q])|\leq{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\frac{c}{n}}. And thus, by the previous result,

We have proved in the beginning of the section that $\frac{1}{T}\operatorname{tr}\Phi$ is bounded and thus we finally conclude that

2.2 Second Equivalent for E​[Q]Edelimited-[]𝑄{\rm E}[Q]

In this section, we show that ${\rm E}[Q]$ can be approximated by the matrix $\bar{Q}$, which we recall is defined as

where $\delta>0$ is the unique positive solution to $\delta=\frac{1}{T}\operatorname{tr}\Phi\bar{Q}$. The fact that $\delta>0$ is well defined is quite standard and has already been proved several times for more elaborate models. Following the ideas of (Hoydis, Couillet and Debbah, 2013), we may for instance use the framework of so-called standard interference functions (Yates, 1995) which claims that, if a map $f:[0,\infty)\to(0,\infty)$, $x\mapsto f(x)$, satisfies $x\geq x^{\prime}\Rightarrow f(x)\geq f(x^{\prime})$, $\forall a>1,af(x)>f(ax)$ and there exists $x_{0}$ such that $x_{0}\geq f(x_{0})$, then $f$ has a unique fixed point (Yates, 1995, Th 2). It is easily shown that $\delta\mapsto\frac{1}{T}\operatorname{tr}\Phi\bar{Q}$ is such a map, so that $\delta$ exists and is unique.

to prove that $|\alpha-\delta|\leq cn^{\varepsilon-\frac{1}{2}}$. To this end, note that, by Cauchy–Schwarz’s inequality,

so that it is sufficient to bound the limsup of both terms under the square root strictly by one. Next, remark that

But at the same time, since $\|(\frac{n}{T}\Phi+\gamma I_{T})^{-1}\|\leq\gamma^{-1}$,

the limsup of which is bounded. We thus conclude that

Similarly, $\alpha$, which is known to be bounded, satisfies

which completes to prove that $|\alpha-\delta|\leq cn^{\varepsilon-\frac{1}{2}}$.

and we have thus proved that $\|{\rm E}[Q]-\bar{Q}\|\leq cn^{-\frac{1}{2}+\varepsilon}$ for some $c>0$.

From this result, along with Corollary 2, we now have that

The evaluation of the second order statistics of the neural network under study requires, beside ${\rm E}[Q]$, to evaluate the more involved form ${\rm E}[QAQ]$, where $A$ is a symmetric matrix either equal to $\Phi$ or of bounded norm (so in particular $\|\bar{Q}A\|$ is bounded). To evaluate this quantity, first write

Of course, since $QAQ$ is symmetric, we may write

which will reveal more practical to handle.

First note that, since $\left\|{\rm E}[Q]-\bar{Q}\right\|\leq Cn^{\varepsilon-\frac{1}{2}}$ and $A$ is such that $\|\bar{Q}A\|$ is bounded, $\|{\rm E}[\bar{Q}AQ]-\bar{Q}A\bar{Q}]\|\leq\|\bar{Q}A\|\|{\rm E}[Q]-\bar{Q}\|\leq C^{\prime}n^{\varepsilon-\frac{1}{2}}$, which provides an estimate for the first expectation. We next evaluate the last right-hand side expectation above. With the same notations as previously, from exchangeability arguments and using $Q=Q_{-}-Q\frac{1}{T}\sigma\sigma^{\sf T}Q_{-}$, observe that

which, reusing $Q=Q_{-}-Q\frac{1}{T}\sigma\sigma^{\sf T}Q_{-}$, is further decomposed as

(where in the previous to last line, we have merely reorganized the terms conveniently) and our interest is in handling $Z_{1}+Z_{1}^{\sf T}+Z_{2}+Z_{2}^{\sf T}+Z_{3}+Z_{3}^{\sf T}+Z_{4}+Z_{4}^{\sf T}$. Let us first treat term $Z_{2}$. Since $\bar{Q}AQ_{-}$ is bounded, by Lemma 4, $\frac{1}{T}\sigma^{\sf T}\bar{Q}AQ_{-}\sigma$ concentrates around $\frac{1}{T}\operatorname{tr}\Phi\bar{Q}AE[Q_{-}]$; but, as $\|\Phi\bar{Q}\|$ is bounded, we also have $|\frac{1}{T}\operatorname{tr}\Phi\bar{Q}AE[Q_{-}]-\frac{1}{T}\operatorname{tr}\Phi\bar{Q}A\bar{Q}|\leq cn^{\varepsilon-\frac{1}{2}}$. We thus deduce, with similar arguments as previously, that

with probability exponentially close to one, in the order of symmetric matrices. Taking expectation and norms on both sides, and conditioning on the aforementioned event and its complementary, we thus have that

with $D=\operatorname{\rm diag}(\{1+\frac{1}{T}\sigma_{i}^{\sf T}Q_{-}\sigma_{i}\})$, the operator norm of which is bounded as O(1). So finally,

We now move to term $Z_{3}+Z_{3}^{\sf T}$. Using the relation $ab^{\sf T}+ba^{\sf T}\preceq aa^{\sf T}+bb^{\sf T}$,

and the symmetrical lower bound (equal to the opposite of the upper bound), where $D_{3}=\operatorname{\rm diag}((\delta-\frac{1}{T}\sigma_{i}^{\sf T}Q_{-i}\sigma_{i})/(1+\frac{1}{T}\sigma_{i}^{\sf T}Q_{-i}\sigma_{i}))$. For the same reasons as above, the first right-hand side term is bounded by $Cn^{\varepsilon-\frac{1}{2}}$. As for the second term, for $A=I_{T}$, it is clearly bounded; for $A=\Phi$, using $\frac{n}{T}\frac{\bar{Q}\Phi}{1+\delta}=I_{T}-\gamma\bar{Q}$, ${\rm E}[Q_{-}A\bar{Q}\Phi\bar{Q}AQ_{-}]$ can be expressed in terms of ${\rm E}[Q_{-}\Phi Q_{-}]$ and ${\rm E}[Q_{-}\bar{Q}^{k}\Phi Q_{-}]$ for $k=1,2$, all of which have been shown to be bounded (at most by $Cn^{\varepsilon}$). We thus conclude that

Finally, term $Z_{4}$ can be handled similarly as term $Z_{2}$ and is shown to be of norm bounded by $Cn^{\varepsilon-\frac{1}{2}}$.

As a consequence of all the above, we thus find that

It is attractive to feel that the sum of the second and third terms above vanishes. This is indeed verified by observing that, for any matrix $B$,

with $D=\operatorname{\rm diag}(1+\frac{1}{T}\sigma_{i}^{\sf T}Q_{-i}\sigma_{i})$, and a similar reasoning is performed to control ${\rm E}[Q_{-}BQ]-{\rm E}[Q_{-}BQ_{-}]$ and ${\rm E}[QBQ_{-}]-{\rm E}[Q_{-}BQ_{-}]$. For $B$ bounded, $\|{\rm E}[Q\frac{1}{T}\Sigma^{\sf T}D\Sigma QBQ]\|$ is bounded as $O(1)$, and thus $\|{\rm E}[QBQ]-{\rm E}[Q_{-}BQ_{-}]\|$ is of order $O(n^{-1})$. So in particular, taking $A$ of bounded norm, we find that

Take now $B=\Phi$. Then, from the relation $AB^{\sf T}+BA^{\sf T}\preceq AA^{\sf T}+BB^{\sf T}$ in the order of symmetric matrices,

The first norm in the parenthesis is bounded by $Cn^{\varepsilon}$ and it thus remains to control the second norm. To this end, similar to the control of ${\rm E}[Q\Phi Q]$, by writing ${\rm E}[Q\Phi Q\Phi Q]={\rm E}[Q\sigma_{1}\sigma_{1}^{\sf T}Q\sigma_{2}\sigma_{2}^{\sf T}Q]$ for $\sigma_{1},\sigma_{2}$ independent vectors with the same law as $\sigma$, and exploiting the exchangeability, we obtain after some calculus that ${\rm E}[Q\Phi Q]$ can be expressed as the sum of terms of the form ${\rm E}[Q_{++}\frac{1}{T}\Sigma_{++}^{\sf T}D\Sigma_{++}Q_{++}]$ or ${\rm E}[Q_{++}\frac{1}{T}\Sigma_{++}^{\sf T}D\Sigma_{++}Q_{++}\frac{1}{T}\Sigma_{++}^{\sf T}D_{2}\Sigma_{++}Q_{++}]$ for $D,D_{2}$ diagonal matrices of norm bounded as $O(1)$, while $\Sigma_{++}$ and $Q_{++}$ are similar as $\Sigma$ and $Q$, only for $n$ replaced by $n+2$. All these terms are bounded as $O(1)$ and we finally obtain that ${\rm E}[Q\Phi Q\Phi Q]$ is bounded and thus

With the additional control on $Q\Phi Q_{-}-Q_{-}\Phi Q_{-}$ and $Q_{-}\Phi Q-Q_{-}\Phi Q_{-}$, together, this implies that {\rm E}[Q\Phi Q]={\rm E}[Q_{-}\Phi Q_{-}]+O_{\|\cdot\|}({\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}n^{-1}}). Hence, for $A=\Phi$, exploiting the fact that $\frac{n}{T}\frac{1}{1+\delta}\Phi\bar{Q}\Phi=\Phi-\gamma\bar{Q}\Phi$, we have the simplification

We have already shown in (11) that $\limsup_{n}\frac{n}{T}\frac{\frac{1}{T}\operatorname{tr}\Phi^{2}\bar{Q}^{2}}{(1+\delta)^{2}}<1$ and thus

which proves immediately Proposition 1 and Theorem 3.

In this section, we evaluate the terms $\Phi_{ab}$ provided in Table 1. The proof for the term corresponding to $\sigma(t)={\rm erf}(t)$ can be already be found in (Williams, 1998, Section 3.1) and is not recalled here. For the other functions $\sigma(\cdot)$, we follow a similar approach as in (Williams, 1998), as detailed next.

The evaluation of $\Phi_{ab}$ for $w\sim\mathcal{N}(0,I_{p})$ requires to estimate

Assume that $a$ and $b$ and not linearly dependent. It is convenient to observe that this integral can be reduced to a two-dimensional integration by considering the basis $e_{1},\ldots,e_{p}$ defined (for instance) by

The case where $a$ and $b$ would be linearly dependent can then be obtained by continuity arguments.

It is worth noticing that this may be more compactly written as

which is minimum for $\angle(a,b)\to-1$ (since $\arccos(-x)\geq 0$ on $$) and takes there the limiting value zero. Hence$\mathcal{I}>0$for$a$and$b$ not linearly dependent.

For $a$ and $b$ linearly dependent, we simply have $\mathcal{I}=0$ for $\angle(a,b)=-1$ and $\mathcal{I}=\frac{1}{2}\|a\|\|b\|$ for $\angle(a,b)=1$.

Since $|t|=\max(t,0)+\max(-t,0)$, we have

Hence, reusing the results above, we have here

Using the identity $\operatorname{acos}(-x)-\operatorname{acos}(x)=2\operatorname{asin}(x)$ provides the expected result.

With the same notations as in the case $\sigma(t)=\max(t,0)$, we have to evaluate