What are the main theorems in Machine (Deep) Learning? Al Rahimi has recently given a very provocative talk in NIPS 2017 comparing current Machine Learning to Alchemy. One of his claims is that we need to get back to theoretical developments, to have simple theorems proving foundational results. 
When he said that, I started looking for the main theorems for ML, but could not find a good reference making sense of the main results. So here is my question: what are the current main mathematical theorems (theory) in ML/DL and what do they prove? I would guess Vapnik's work would go somewhere here. As an extra, what are the main theoretical open problems?
 A: As I wrote in the comments, this question seems too broad to me, but I'll make an attempt to an answer. In order to set some boundaries, I will start with a little math which underlies most of ML, and then concentrate on recent results for DL.

The bias-variance tradeoff is referred to in countless books, courses, MOOCs, blogs, tweets, etc. on ML, so we can't start without mentioning it:
$$\mathbb{E}[(Y-\hat{f}(X))^2|X=x_0]=\sigma_{\epsilon}^2+\left(\mathbb{E}\hat{f}(x_0)-f(x_0)\right)^2+\mathbb{E}\left[\left(\hat{f}(x_0)-\mathbb{E}\hat{f}(x_0)\right)^2\right]=\text{Irreducible error + Bias}^2 \text{ + Variance}$$
Proof here: https://web.stanford.edu/~hastie/ElemStatLearn/

The Gauss-Markov Theorem (yes, linear regression will remain an important part of Machine Learning, no matter what: deal with it) clarifies that, when the linear model is true and some assumptions on the error term are valid, OLS has the minimum mean squared error (which in the above expression is just $\text{Bias}^2 \text{ + Variance}$) only among the unbiased linear estimators of the linear model. Thus there could well be linear estimators with bias (or nonlinear estimators) which have a better mean square error, and thus a better expected prediction error, than OLS. And this paves the way to all the regularization arsenal (ridge regression, LASSO, weight decay, etc.) which is a workhorse of ML. A proof is given here (and in countless other books):
https://www.amazon.com/Linear-Statistical-Models-James-Stapleton/dp/0470231467
Probably more relevant to the explosion of regularization approaches, as noted by Carlos Cinelli in the comments, and definitely more fun to learn about, is the James-Stein theorem. Consider $n$ independent, same variance but not same mean Gaussian random variables:
$$X_i|\mu_i\sim \mathcal{N}(\theta_i,\sigma^2), \quad i=1,\dots,n$$
in other words, we have an $n-$components Gaussian random vector $\mathbf{X}\sim \mathcal{N}(\boldsymbol{\theta},\sigma^2I)$. We have one sample $\mathbf{x}$ from $\mathbf{X}$ and we want to estimate $\boldsymbol{\theta}$. The MLE (and also UMVUE) estimator is obviously $\hat{\boldsymbol{\theta}}_{MLE}=\mathbf{x}$. Consider the James-Stein estimator
$$\hat{\boldsymbol{\theta}}_{JS}= \left(1-\frac{(n-2)\sigma^2}{||\mathbf{x}||^2}\right)\mathbf{x} $$
Clearly, if $(n-2)\sigma^2\leq||\mathbf{x}||^2$, $\hat{\boldsymbol{\theta}}_{JS}$ shrinks the MLE estimate towards zero. The James-Stein theorem states that for $n\geq4$, $\hat{\boldsymbol{\theta}}_{JS}$ strictly dominates $\hat{\boldsymbol{\theta}}_{MLE}$, i.e., it has lower MSE $\forall \ \boldsymbol{\theta}$. Pheraps surprisingly, even if we shrink towards any other constant $\boldsymbol{c}\neq \mathbf{0}$, $\hat{\boldsymbol{\theta}}_{JS}$ still dominates $\hat{\boldsymbol{\theta}}_{MLE}$. Since the $X_i$ are independent, it may seem weird that, when trying to estimate  the height of three unrelated persons, including a sample from the number of apples produced in Spain, may improve our estimate on average. The key point here is "on average": the mean square error for the simultaneous estimation of all the components of the parameter vector is smaller, but the square error for one or more components  may well be larger, and indeed it often is, when you have "extreme" observations.
Finding out that MLE, which was indeed the "optimal" estimator for the univariate estimation case, was dethroned for multivariate estimation, was quite a shock at the time, and led to a great interest in shrinkage, better known as regularization in ML parlance. One could note some similarities with mixed models and the concept of "borrowing strength": there is indeed some connection, as discussed here
Unified view on shrinkage: what is the relation (if any) between Stein's paradox, ridge regression, and random effects in mixed models?
Reference: James, W., Stein, C., Estimation with Quadratic Loss. Proceedings of the Fourth Berkeley Symposium on Mathematical Statistics and Probability, Volume 1: Contributions to the Theory of Statistics, 361--379, University of California Press, Berkeley, Calif., 1961

Principal Component Analysis is key to the important topic of dimension reduction, and it's based on the Singular Value Decomposition: for each $N\times p$ real matrix $X$ (although the theorem easily generalizes to complex matrices) we can write
$$X=UDV^T$$
where $U$ of size $N \times p$ is orthogonal, $D$ is a $p \times p$ diagonal matrix with nonnegative diagonal elements and $U$ of size $p \times p$ is again orthogonal. For proofs and algorithms on how to compute it see: Golub, G., and Van Loan, C. (1983), Matrix computations, John Hopkins University press, Baltimore.

Mercer's theorem is the founding stone for a lot of different ML methods: thin plate splines, support vector machines, the Kriging estimate of a Gaussian random process, etc. Basically, is one of the two theorems behind the so-called kernel trick. Let $K(x,y):[a,b]\times[a,b]\to\mathbb{R}$ be a symmmetric continuous function or kernel. if $K$ is positive semidefinite, then it admits an orthornormal basis of eigenfunctions corresponding to nonnegative eigenvalues:
$$K(x,y)=\sum_{i=1}^\infty\gamma_i \phi_i(x)\phi_i(y)$$
The importance of this theorem for ML theory is testified by the number of references it gets in famous texts, such as for example Rasmussen & Williams text on Gaussian processes.
Reference: J. Mercer, Functions of positive and negative type, and their connection with the theory of integral equations. Philosophical Transactions of the Royal Society of London. Series A, Containing Papers of a Mathematical or Physical Character, 209:415-446, 1909
There is also a simpler presentation in Konrad Jörgens, Linear integral operators, Pitman, Boston, 1982.

The other theorem which, together with Mercer's theorem, lays out the theoretical foundation of the kernel trick, is the representer theorem. Suppose you have a sample space $\mathcal{X}$ and a symmetric positive semidefinite kernel $K: \mathcal{X} \times \mathcal{X}\to \mathbb{R}$. Also let $\mathcal{H}_K$ be the RKHS associated with $K$. Finally, let $S=\{\mathbb{x}_i,y_i\}_{i=1}^n$ be a training sample. The theorem says that among all functions $f\in \mathcal{H}_K$, which all admit an infinite representation in terms of eigenfunctions of $K$ because of Mercer's theorem, the one that minimizes the regularized risk always has a finite representation in the basis formed by the kernel evaluated at the $n$ training points, i.e.
$$\min_{f \in \mathcal{H}_K} \sum_{i=1}^n L(y_i,f(x_i))+\lambda||f||^2_{\mathcal{H}_K}=\min_{\{c_j\}_1^\infty} \sum_{i=1}^n L(y_i,\sum_j^\infty c_j\phi_j(x_i))+\lambda\sum_j^\infty \frac{c_j^2}{\gamma_j}=\sum_{i=1}^n\alpha_i K(x,x_i)$$
(the theorem is the last equality). References: Wahba, G. 1990, Spline Models for Observational Data, SIAM, Philadelphia.

The universal approximation theorem has been already cited by user Tobias Windisch and is much less relevant to Machine Learning than it is to functional analysis, even if it may not seem so at a first glance. The problem is that the theorem only says that such a network exists, but:

*

*it doesn't give any correlation between the size $N$ of the hidden layer and some measure of complexity of the target function $f(x)$, such as for example Total Variation. If $f(x)=\sin(\omega x):[0,2\pi]\to[-1,1]$ and the $N$ required for a fixed error $\epsilon$ growed exponentially with $\omega$, then single hidden layer neural networks would be useless.

*it doesn't say if the network $F(x)$ is learnable. In other words assume that given $f$ and $\epsilon$, we know that a size $N$ NN will approximate $f$ with the required tolerance in the hypercube. Then by using training sets of size $M$  and a learning procedure such as for example back-prop, do we have any guarantee that by increasing $M$ we can recover $F$?

*finally, and worse of them all, it doesn't say anything about the prediction error of neural networks. What we're really interested in is an estimate of the prediction error, at least averaged over all training sets of size $M$. The theorem doesn't help in this respect.

A smaller pain point with the Hornik's version of this theorem is that it doesn't hold for ReLU activation functions. However, Bartlett has since proved  an extended version which covers this gap.

Until now, I guess all the theorems I considered were well-known to anybody. So now it's time for the fun stuff :-) Let's see a few Deep Learning theorems:
Assumptions:

*

*the deep neural network $\Phi(X,W)$ (for fixed $W$, $\Phi_W(X)$ is the function which associates the inputs of the neural network with its outputs) and the regularization loss $\Theta(W)$ are both sums of positively homogeneous functions of the same degree

*the loss function $L(Y,\Phi(X,W))$ is convex and at least once differentiable in $X$, in a compact set $S$
Then:

*

*any local minimum for $L(Y,\Phi(X,W))+\lambda\Theta(W)$ such that a subnetwork of $\Phi(X,W)$ has zero weights, is a global minimum (Theorem 1)

*above a critical network size, local descent will always converge to a global minimum from any initialization (Theorem 2).

This is very interesting: CNNs made only of convolutional layers, ReLU, max-pooling, fully connected ReLU and linear layers are positively homogenous functions, while if we include sigmoid activation functions, this isn't true anymore, which may partly explain the superior performance in some applications of ReLU + max pooling with respect to sigmoids. What's more, the theorems only hold if also $\Theta$ is positively homogeneous in $W$ of the same degree as $\Phi$. Now, the fun fact is that  $l_1$ or $l_2$ regularization, although positively homogeneous, don't have the same degree of $\Phi$ (the degree of $\Phi$, in the simple CNN case mentioned before, increases with the number of layers). Instead, more modern regularization methods such as batch normalization and path-SGD do correspond to a positively homogeneous regularization function of the same degree as $\Phi$, and dropout, while not fitting this framework exactly, holds strong similarities to it. This may explain why, in order to get high accuracy with CNNs, $l_1$ and $l_2$ regularization are not enough, but we need to employ all kinds of devilish tricks, such as dropout and batch normalization! To the best of my knowledge, this is the closest thing to an explanation of the efficacy of batch normalization, which is otherwise very obscure, as correctly noted by Al Rahimi in his talk.
Another observation that some people make, based on Theorem 1, is that it could explain why ReLU work well, even with the problem of dead neurons. According to this intuition, the fact that, during training, some ReLU neurons "die" (go to zero activation and then never recover from that, since for $x<0$ the gradient of ReLU is zero) is "a feature, not a bug", because if we have reached a minimum and a full subnetwork has died, then we're provably reached a global minimum (under the hypotheses of Theorem 1). I may be missing something, but I think this interpretation is far-fetched. First of all, during training ReLUs can "die" well before we have reached a local minimun. Secondly, it has to be proved that when ReLU units "die", they always do it over a full subnetwork: the only case where this is trivially true is when you have just one hidden layer, in which case of course each single neuron is a subnetwork. But in general I would be very cautious in seeing "dead neurons" as a good thing.
References:
B. Haeffele and R. Vidal, Global optimality in neural network training,
In IEEE Conference on Computer Vision and Pattern Recognition,
2017.
B. Haeffele and R. Vidal. Global optimality in tensor factorization,
deep learning, and beyond, arXiv, abs/1506.07540, 2015.

Image classification requires learning representations which are invariant (or at least robust, i.e., very weakly sensitive) to various transformations such as  location, pose, viewpoint, lighting, expression, etc. which are commonly present in natural images, but do not contain info for the classification task. Same thing for speech recognition: changes in pitch, volume, pace, accent. etc. should not lead to a change in the classification of the word. Operations such as convolution, max pooling, average pooling, etc., used in CNNs, have exactly this goal, so intuitively we expect that they would work for these applications. But do we have theorems to support this intuition? There is a vertical translation invariance theorem, which, notwithstanding the name, has nothing to do with translation in the vertical direction, but it's basically a result which says that features learnt in following layers get more and more invariant, as the number of layers grows. This is opposed to an older horizontal translation invariance theorem which however holds for scattering networks, but not for CNNs.
The theorem is very technical, however:

*

*assume $f$ (your input image) is square-integrable

*assume your filter commutes with the translation operator $T_t$, which maps the input image $f$ to a translated copy of itself $T_t f$. A learned convolution kernel (filter) satisfies this hypothesis.

*assume all filters, nonlinearities and pooling in your network satisfy a so-called weak admissibility condition, which is basically some sort of weak regularity and boundedness conditions. These conditions are satisfied by learned convolution kernel (as long as some normalization operation is performed on each layer), ReLU, sigmoid, tanh, etc, nonlinearities, and by average pooling, but not by max-pooling. So it covers some (not all) real world CNN architectures.

*Assume finally that each layer $n$ has a pooling factor $S_n> 1$, i.e., pooling is applied in each layer and effectively discards information. The condition $S_n\geq 1 $ would also suffice for a weaker version of the theorem.

Indicate with $\Phi^n(f)$ the output of layer $n$ of the CNN, when the input is $f$. Then finally:
$$\lim_{n\to\infty}|||\Phi^n(T_f f)-\Phi^n(f)|||=0$$
(the triple bars are not an error) which basically means that each layer learns features which become more and more invariant, and in the limit of an infinitely deep network we have a perfectly invariant architecture. Since CNNs have a finite number of layers, they're not perfectly translation-invariant, which is something well-known to practitioners.
Reference: T. Wiatowski and H. Bolcskei, A Mathematical Theory of Deep Convolutional Neural Networks for Feature Extraction, arXiv:1512.06293v3.

To conclude, numerous bounds for the generalization error of a Deep Neural Network based on its Vapnik-Chervonkensis dimension or on the Rademacher complexity grow with the number of parameters (some even exponentially), which means they can't explain why DNNs work so well in practice even when the number of parameters is considerably larger than the number of training samples. As a matter of fact, VC theory is not very useful in Deep Learning.
Conversely, some results from last year bound the generalization error of a DNN classifier with a quantity which is independent of the neural network's depth and size, but depends only on the structure of the training set and the input space. Under some pretty technical assumptions on the learning procedure, and on the training set and input space, but with very little assumptions on the DNN (in particular, CNNs are fully covered), then with probability at least $1-\delta$, we have
$$\text{GE} \leq \sqrt{2\log{2}N_y\frac{\mathcal{N_{\gamma}}}{m}}+\sqrt{\frac{2\log{(1/\delta)}}{m}}$$
where:

*

*$\text{GE}$ is the generalization error, defined as the difference between the expected loss (the average loss of the learned classifier on all possible test points) and the empirical loss (just the good ol' training set error)

*$N_y$ is the number of classes

*$m$ is the size of the training set

*$\mathcal{N_{\gamma}}$ is the covering number of the data, a quantity related to the structure of the input space and to the the minimal separation among points of different classes in the training set. Reference:

J. Sokolic, R. Giryes, G. Sapiro, and M. Rodrigues. Generalization
error of invariant classifiers. In AISTATS, 2017
A: The Kernel Trick is a general idea that's used in a lot of places, and comes from a lot of abstract maths about Hilbert Spaces. Way too much theory for me to type (copy...) out into an answer here, but if you skim through this you can get a good idea of its rigorous underpinnings:
http://www.stats.ox.ac.uk/~sejdinov/teaching/atml14/Theory_2014.pdf
A: My favourite one is the Kraft inequality.

Theorem: For any description method $C$ for finite alphabet $A = \{1,\dots, m\}$, the code word lengths $L_C(1), \dots, L_C(2)$ must satisfy the inequality $\sum_{x \in A} 2 ^{-L_C(x)} \leq 1$.

This inequality relates compression with probability densities: given a code, the length of an outcome represented by that code is the negative log probability of a model identified by the code.
Further, the no free lunch theorem for machine learning has a less well known sibling the no hyper compression theorem, which states that not all sequences can be compressed.
A: I wouldn't call it a main theorem, but I think the following (sometimes referred to as the Universal approximation theorem) is an interesting (and at least for me surprising) one as it states the approximative power of feed-forward neural networks. 
Theorem:
Let $\sigma$ be a nonconstant and monotinically-increasing continous function. For any continuos function $f:[0,1]^m\to\mathbb{R}$ and any $\epsilon>0$, there exist an integern $N$ and an multilayer perceptron $F$ with one hidden layer having $N$ Neurons that has $\sigma$ as activation function so that
$$|F(x)-f(x)|\le\epsilon$$
for all $x\in[0,1]^m$.
Of course, as this is a statement on existence, its impact for practitioners is negligible. 
A proof can be found in  Hornik, Approximation Capabilities of Muitilayer Feedforward Networks, Neural Networks 4 (2), 1991,
A: A nice post focusing on this question (specifically deep learning rather than general machine learning theorems) is here:
https://medium.com/mlreview/modern-theory-of-deep-learning-why-does-it-works-so-well-9ee1f7fb2808
It gives an accessible summary of the main emerging theorems for the ability of deep neural networks to generalize so well. 
A: I think the following theorem that you allude to is considered to be pretty fundamental in statistical learning.
Theorem (Vapnik and Chervonenkis, 1971) Let $H$ be a
hypothesis class of functions from a domain $X$ to $\{0, 1\}$ and let the loss function
be the $0 − 1$ loss. Then, the following are equivalent:


*

*$H$ has the uniform convergence property.

*$H$ is PAC learnable.

*$H$ has a finite VC-dimension.


Proved in a quantitive version here: 
V.N. Vapnik and A.Y. Chervonenkis: On the uniform convergence of relative
frequencies of events to their probabilities. Theory of Probability and its Applications,
16(2): 264–280, 1971.
The version formulated above along with an nice exposition of other results from learning theory is available here:
Shalev-Shwartz, Shai, and Shai Ben-David. Understanding machine learning: From theory to algorithms. Cambridge university press, 2014.
