Intelligent Matrix Exponentiation
Authors: Thomas Fischbacher
Iulia Maria Comsa
Krzysztof Potempa
Moritz Firsching
Luca Versari
Jyrki Alakuijala
Venue:  (2021)
Keywords: matrix exponential
tensor methods
supervised learning
domain extrapolation
certified robustness

Abstract: We present a novel machine learning architecture that uses a single high-dimensional nonlinearity consisting of the exponential of a single input-dependent matrix. The mathematical simplicity of this architecture allows a detailed analysis of its behaviour, providing robustness guarantees via Lipschitz bounds. Despite its simplicity, a single matrix exponential layer already provides universal approximation properties and can learn and extrapolate fundamental functions of the input, such as periodic structure or geometric invariants. This architecture outperforms other general-purpose architectures on benchmark problems, including CIFAR-10, using fewer parameters.

Under review as a conference paper at ICLR 2021
INTELLIGENT MATRIX EXPONENTIATION
Anonymous authors Paper under double-blind review
ABSTRACT
We present a novel machine learning architecture that uses a single high- dimensional nonlinearity consisting of the exponential a single input- matrix. The mathematical simplicity of this architecture allows a detailed analy- sis of its hehaviour, providing robustness guarantees via Lipschitz bounds. De- spite its simplicity, single matrix exponential layer already provides universal approximation properties and can learn and extrapolate fundamental functions of the input, such as periodic structure or geometric invariants. This architecture out- performs other general-purposs architectures on benchmark problems, including CIFAR-10 10, using fewer parameters.
INTRODUCTION
Deep neural networks (DNNs) synthesize highly complex functions by composing a large num- ber of neuronal units, each featuring basic and usually 1-dimensional nonlinear activation func- tion f IR_ + R While highly successful in practice, this approach also has disadvantages. In conventional DNN, any activations only ever get combined through summation. This means that such network requires an increasing number of parameters to approximate more complex functions even simple as multiplication. Parameter-wist this approach of composing simple functions does not scale efficiently.
An alternative to the composition of many 1-dimensional functions is using a simple higher- dimensional nonlinear function Im + RR. A single multidimensional nonlinearity may be desirable because it could express more complex relationships between input features with poten- tially fewer parameters and fewer mathematical operations.
The matrix exponential stands out as a promising but overlooked candidate for a higher-dimensionall nonlinearity for machine learning models. The matrix exponential is smooth function that appears in the solution to one of the simplest differential equations that can yield desirable mathematical properties: d/dty((t = My((), with the solution y(t) exp( exp(M))y(OO), where M is constant matrix. The matrix exponential plays a prominent role in the theory of Lie groups, an algebraic structure widely used throughout many branches of mathematics and science.
unique advantage of the matrix exponential its natural ability to represent oscillations and exponential decay, which becomes apparent if we decompose the matrix to be exponentiated into a symmetric and an antisymmetric component. The exponential of an antisymmetric matrix, whose nonzero eigenvalues are always imaginary, generates superposition of periodic oscillations, whereas the exponential of symmetric matrix, which has real eigenvalues, expresses an exponential growth or decay. This fact, especially the ability represent periodic functions, gives the M-layer the possibility to extrapolate beyond its training domain. Many real-world phenomena contain some degree of periodicity and can therefore benefit from this feature. contrast, the functions typically used in conventional DNNs approximate the target function locally and are therefore unable to scale well outside the boundaries of the training data for such problems.
Based on these insights, we propose a novel architecture for supervised learning whose core element a single layer (henceforth referred to as "M-layer"), that computes a single matrix exponential, where the matrix to be exponentiated is affine function of the input features. We show that the M-layer has universal approximator properties and allows closed-form. per-example bounds for ro- bustness. We demonstrate the ability of this architecture learn multivariate polynomials, such as matrix determinants, and to generalize periodic functions beyond the domain of the input without any feature engineering. Furthermore, the M-layer achieves results comparable eccently-roposee

---

Under review as a conference paper at ICLR 2021
non-specialized architectures image recognition datasets. We provide TensorFlow code that im- plements the M-layer in the Supplementary Material.
RELATED WORK
Neuronal units with more complex activation functions have been proposed. One such example are sigma-pi units (Rumelhart et al., 1986), whose activation function is the weighted sum of products of its inputs. More recently, neural arithmetic logic units have been introduced (Trask et al., 2018), which can combine inputs using multiple arithmetic operators and generalize outside the domain of the training data. In contrast with these architectures, the M-layer is not based on neuronal units with multiple inputs, but uses single matrix exponential as its nonlinear mapping function. Through the matrix exponential, the M-layer can casily learn mathematical operations more complex than addition, but with simpler architecture. In fact, as shown in Section 3.3, the M-layer can be regarded as a generalized sigma-pi network with built-in architecture search, in the sense that learns by itself which arithmetic graph should be used for the computation.
Architectures with higher-dimensional nonlinearities are also already used. The softmax function an example of widely-used such nonlinear activation function that, like the M-layer, has extra mathematical structure. For example, a permutation of the softmax imputs produces corresponding permutation of the outputs. Maxout networks also act on multiple units and have been successful in combination with dropout (Goodfellow et al., 2013). In radial basis networks (Park & Sandberg, 1991), each hidden unit computes nonlinear function of the distance between its own learned cen- troid and a single point represented by vector of input coordinates. Capsule networks (Sabour et al., 2017) are another recent example of multidimensional nonlinearities. Similarly, the M-layer uses the matrix exponential as single high-dimensional nonlinearity. This mapping satisfies non- trivial mathematical identities, for some of which we here elucidate their application potential in machine learning context. Feature crosses of second and higher order have also been explored on top of convolutional layers ((Lin et al., 2015; Lin & Maji, 2017; Li et al., 2017; Engin et al., 2018; Koniusz al., 2018). In this line of research, (Li al., 2017) introduces the use of matrix operation (cither matrix logarithm or matrix square root) on top of a convolutional layer with higher order feature crosses. The M-Layer differs from this kind of approach in two ways: first, uses ma- trix exponentiation to produce feature crosses, and not to improve trainability on top of other ones; secondly, computes a matrix operation of an arbitrary matrix (not necessarily symmetrie positive semidefinite), which allows more expressibility, as seen for example in Section 3.3.
The main differences between the M-layer and other architectures that can utilize feature-crosser are the M-layer's ability to also model non-polynomial dependencies (such as a cosine) and the built- in competition-for wow.ino---oo-otll--ooooooooo "circuit breadboard" analogy.
regularization, which we will explain in Section 3.3 through a
Matrix exponentiation has a natural alternative interpretation in terms of an ordinary differential equation (ODE). As such, the M-layer can be compared to other novel ODE-related architectures, such as neural ordinary diffcrential equations (NODE) (Chen et al., 2018) and their augmented extensions (ANODE) (Dupont et al., 2019). We discuss this in Section 3.6.
Existing approaches to certifying the robustness of neural networks can be split into two different categories. Some approaches (Peck et al., 2017) mathematically analyze network layer by layer, providing bounds on the robustness of each layer, that then get multiplied together. This kind of approach tends to give fairly loose bounds, due to the inherent tightness loss from composing upper bounds. Other approaches (Singh et al., 2018; 2019) use abstract interpretation on the evaluation of the network provide empirical robustness bounds. In contrast, using the fact that the M-layer architecture has a single layer, in Section 3.7 we obtain direct bound on the robustness on the whole network by analyzing the explicit formulation of the computation.
3 ARCHITECTURE
We start this section by refreshing the definition of the matrix exponential. We then define the pro- posed M-layer architecture and explain its ability to learn particular functions such as polynomials and periodic functions. Finally, provide closed-form per-example robustness guarantees.
2

---

Under review as a conference paper at ICLR 2021
3.1 MATRIX EXPONENTIATION
The exponential of a square matrix M is defined as: exp(M) - ÆT MM. The matrix power M is defined inductively as M = I, M = M MM,, using the associativity of the matrix product; is not an element-wise matrix operation. Note that the expansion of exp(M) (M) in Eq. (3.1) is finite for nilpotent matrices. A matrix M is called nilpotent if there exists positive integer & such that M O. Strictly upper triangular matrices are a canonical example.
Multiple algorithms for computing the matrix exponential efficiently have been proposed (Moler & Van Loan, 2003). We provide an implementation of the M-layer that employs either the TensorFlow function tf.. 1 linalg.. expm, which uses the scaling and squaring method combined with the Padé approximation (Higham, 2005), or a simpler and faster approximation of the matrix exponential that
only squares a fixed number of times, using the formula exp(M) = (I (I+ 2 ) This approximation follows from the equivalent definition of cXp(M) lim, [ )", and has (knn--) time complexity, where n is the total number of elements M.
M
3.2 M-LAYER DEFINITION
At the core of the proposed architecture (Figure 1) is an M-layer that computes a single matrix exponential, where the matrix to be exponentiated is an affine function of all of the input features. In other words, an M-layer replaces an entire stack of hidden layers DNN.
We exemplify the architecture as applied to standard im- age recognition dataset, but this formulation is applicable to any other type of problem by adapting the relevant in- put indices. In the following equations, generalized Ein- stein summation is performed over all right- side in- dices not seen on the left- left-hand side.
exp(-)
Consider example input image, encoded 3-index array Xyze. where y. x and c are the row index, column index and color channel index, respectively. The matrix M to be exponentiated is obtained as follows, using the trainable parameters Tajts. U Uazye and Bik:
Uz TUw exp(TUx)
Sexp(TUx)
Figure 1: M-layer architecture diagram.
Mfk Byx Pagh UayzeX you
(1)
First, X is projected linearly to d-dimensional latent feature embedding space by Üazze- Then, the 3-index tensor Tajk maps each such latent feature to an n X n matrix. Finally, a bias matrix Byx is added to the peature-weighted sum of matrices. The result is a matrix indexed by row and column indices j and k. In order to allow for exponentiation, M must be a square matrix. is possible to contract the tensors T and Ü in order to simplify the architecture formula, but partial tensor factorization provides regularization by reducing the parameter count. An output pm is obtained as follows, using the trainable parameters Smjj and Vm:
Pm Vn + Smyk exp(MM
The matrix exp(M), indexed by row and column indices and k: in the same way M, is projected linearly by the 3-index tensor Smijk to obtain a h--dimensional output vector. The bias-vector V... turns this linear mapping into an affine mapping. The resulting vector may be interpreted as accu- mulated per-class evidence and may then be mapped to a vector of probabilities via softmax.
Training can thus be done conventionally. by minimizing a loss function such as the L2 norm or the cross-entropy with softmax, using backpropagation through matrix exponentiation. In general, there is no simple analytical expression for the derivative of cxp(M), but one can instead com- pute derivatives for each step of the algorithm chosen to compute the exponential. For example, f.. 1 i naa 1 g. expm uses the scaling-and-squu steps are all differentiable. We discuss the topic in more detail Appendix A.2.
algorithm with a Padé approximation, whose
The nonlinearity of the M-layer architecture is provided by the Rd + R² mapping + V. + Smnk exp( (M) jkc The count of trainable parameters of this component dn? n² n²h + h. This
3

---

Under review as a conference paper at ICLR 2021
count comes from summing the dimensions of Tajk, Byk: Smiks and Vn respectively. We note that this architecture has some redundancy in its parameters, as one can freely multiply the T and U tensors by a d x d real matrix and, respectively, its inverse, while preserving the computed function. Similarly, it is possible to multiply each of the x n parts of the tensors T and S, as well as B, by both an n x n matrix and its inverse. In other words, any pair of real invertible matrices of sizes d x and X n can be used produce a new parametrization that still computes the same function.
3.3 FEATURE CROSSES AND UNIVERSAL APPROXIMATION
A key property of the M-layer is its ability to generate arbitrary exponential-polyyoniin combina- tions of the input features. For classification problems, M-layer architectures are a superset of mul- tivariate polynomial classifiers, where the matrix size constrains the complexity of the polynomial while at the same time not uniformly constraining its degree. other words, simple multivariate polynomials of high degree compete against complex multivariate polynomials of low degree.
We provide universal approximator proof for the M-layer in Appendix A.I, which relies on its ability to express any multivariate polynomial in the input features if sufficiently large matrix size is used. We provide here an example that illustrates how feature crosses can be generated through the matrix exponential. Consider dataset with the feature vector (øu, 01, o2) given by the Ü æ tensor contraction, where the relevant quantities for the final classification of an example 60. 01. 22, Φ1Φ1. and Φ1 01 To learn this dataset, we look for an exponentiated matrix that makes precisely these quantities available to be weighted by the trainable tensor S. To do this, we define three 7 matrices Tojk. Tikk, and T2jk. as To01 Tio2 203- 1, To24 T225 T256 3, and 0 otherwise. We then define the matrix M as:
assumed to be
6n di do
0
U
U
0
0
0
200
202
M To + Ti T2
0
0
0
0
302
0
Note that M is nilpotent, as M = 0. Therefore, we obtain the following matrix exponential, which contains the desired quantities in its leading row:
du
62
0
0
0
0
0
1 M2 + M3 exp(M) I + M
20p
200
3e
2
0
0
0
1
0
1
0
0
0
0
0
0
300
The same technique can be employed to encode any polynomial in the input features using a n n matrix, where is one unit larger than the total number of features plus the intermediate and final products that need to be computed. The matrix size can be seen as regulating the total capacity of the model for computing different feature crosses. With this intuition, "cireuit breadboard" for wiring up arbitrary polynomials.
read the matrix
3.4 FEATURE PERIODICITY
While the M-layer can express a wide range of functions using the exponential of nilpotent matrices, non-nilpotent matrices also have important applications, such as learning the periodicity of input features. This is a problem where conventional DNNs struggle, as they cannot naturally generalize beyond the distribution of the training data. To illustrate how matrix exponentials can naturally fit periodic dependency on input features, without requiring an explicit specification of the periodic nature of the data, consider the matrix M, (2 (2-)). We have exp(tM,) sin wet which is a 2d rotation by an angle of wl and thus periodic in / with period 2 2n/w. This setup can fit functions that have an arbitrary period. Moreover, this representation of periodicity naturally extrapolates well when going beyond the range of the initial numerical data.
wt -sin wt we ).
3.5 CONNECTION TO LIE GROUPS
The M-layer has a natural connection to Lie groups. Lie groups can be thought of as a model of continuous symmetries of system such as rotations. There is large body of mathematical theory
4

---

Under review as a conference paper at ICLR 2021
and tools available to study the structure and properties of Lie groups (Gilmore, 2008; 2012), which may ultimately also help for model interpretability.
Every Lie group has associated Lie algebra, which can be understood as the space of the small perturbations with which it is possible to generate the elements of the Lie group. As an example, the set of rotations of 3-dimensional space forms Lie group; the corresponding algebra can be understood as the set of rotation axes in 3 dimensions. Lie groups and algebras can be represented using matrices, and by computing matrix exponential one can map elements of the algebra to elements of the group.
In the M-layer architecture, the role of the 3-index tensor / is to form a matrix whose entries are affine functions of the input features. The matrices that compose T can be thought of as generators of Lie algebra. Building M corresponds selecting a Lie algebra element. Matrix exponentiation then computes the corresponding Lie group element. As rotations are periodic and one of the sim- plest forms of continuous symmetries, this perspective is useful for understanding the ability of the M-layer to learn periodicity in input features.
The connection Lie groups is important because can provide techniques for the regularization of the M-layer. It is, for example, possible to impose constraint on the maximal depth of feature crosses by taking generator matrix tensor T that itself is tensor product with a nilpotent matrix algebra. Other Lie-group- dic-group-inspired regularization techniques will be discussed in future work.
3.6 DYNAMICAL SYSTEMS INTERPRETATION
Recent work has proposed dynamical systems interpretation of some DNN architectures. The ODE that is provided
NODE architecture (Chen et al., 2018) uses a nonlinear and not by trainable neural units, and computes the time evolution of vector that is constructed from the input features. This section discusses a similar interpretation of the M-layer.
Consider an M-layer with Í defined as To12 = Ti20 = T201 = +1, Teno = Tiez = Tozi = 11 and 0 otherwise, with Ü as the 3 x 3 identity matrix, and with B = O. Given an input vector a, Plugging M into the linear and time invariant eo
the corresponding matrix M is then (-a2 -ap 0
-a:
(LTI) ODE d/ dt (t) = MY (t). we can observe that the ODE describes a rotation around the axis defined by Moreover, solution this ODE is given by (tt) (Xp(tM Thus, by choosing Smik- = Y()) m j and 0 otherwise, the above M-layer can be understood as applying a rotation with input dependent angular velocity some basis vector over a unit time interval.
More generally, we can consider the input features to provide affine parameters that define a time- invariant linear ODE, and the output of the M-layer be affine function of vector that has evolved under the ODE over a unit time interval. In contrast, the NODE architecture uses a non- linear ODE that is not input dependent, which gets applied to input-dependent feature vector.
3.7 CERTIFIED ROBUSTNESS
We show that the M-layer allows novel proof technique to produce closed-form expressions for guaranteed robustness bounds. For any matrix norm [|-1|, we have (Roger & Charles, 1994):
||exp(X + +Y exp( (X)l s [|Y]] exp([I)]) exp([(x]))
We also make use of the fact that || < V M la for any n x n matrix, where ||-|| is the Frobenius norm and I| the 2-norm of matrix. We recall that the Frobenius norm of matrix is equivalent to the 2-norm of the vector formed from the matrix entries.
Let M be the matrix be exponentiated corresponding to a given input example x, and let M' be the deviation to this matrix that corresponds to an input deviation of i. i.e. M + M' is the matrix corresponding to input example x + x. Given that the mapping between and M is linear, there is mor-model constant d,,, such that ||M''] I 6, l The 2-norm of the difference between the outputs can bound as follows:
| la ||||||||expMM + M¹) - exp(Ai) S / + M' exp(M))l
S||2 M |lx[i) (1)) (xx()111 S V || exp( |=]1|=0 exp(
5

---

Under review as a conference paper at ICLR 2021
10
10
-5-
110
-10
(a) RBF-SVM
(b) ReLU DNN
(c) M-layer
Figure? 2: Classification boundaries, training and test sets on a "Swiss roll" classification task.
where ||S||2 is computed by considering S a h x n rectangular matrix, and the first inequality follows from the fact that the tensor multiplication by S can be considered matrix-vector multipli- cation between S and the result of matrix exponential seen as a vector.
This inequality allows to compute the minimal T.x change required in the input given the difference between the amount of accumulated evidence between the most likely class and other classes. More- over, considering that || is bounded from above, for example by I in the case of CIFAR-10, we can obtain Lipschitz bound by replacing the exp(6mn |=][)>)>)))))))))))))))))))))| term with a exp(ô..) term.
alternative way of studying robustness for specific example involves the diagonalization of the matrix that corresponds that example. In particular, (BCH) formula to understand the neighborhood of the example. For example, if we regularize the tensor T to punish non-commutativity of generator matrices, expect the terms in the BCH for- mula to decrease in magnitude according to the number of commutators, so we can obtain an under- standing of the non-linear effects in the neighborhood of each specific example. To our knowledge, no equivalent type of analysis exists for ReLU networks.
use the
RESULTS
In this section, we show the performance of the M-layer on multiple benchmarks, in comparison chosen based on
with more traditional architectures. For each experiment, the hyperparameters performance. The DNNs use uniform Glorot initialization (Glorot & Bengio, 2010). The M-layers initialized using normal distributions (u 0) for the spiral, periodicity and determinants exper- iments (a {1, 0.01, 1}. respectively), and uniformly in [-0.05, 0.05) for the image experiments. regularization term is used on the M-layer for the spiral (L2 10-3) )) determinants (L2 10-1) and image (L1 10-3) experiments. The selected number of epochs, batch size and learning rate are: 1000, 16, 10-2 (spiral); 1000, 128/16 16 (M-layer/RellJU), 10-5/-4 (periodicity): 256, 32, 10-3 (determinants); 150, 32, 10-3 (images). The learning rate is reduced on plateaus and early stopping is used in all experiments except periodicity. RMSprop optimizers are used for all experiments ex- cept for images, where SGD with momentum 0.9 is used. The loss function is the cross-entropy with softmax (spiral, images) or the mean squared error (MSE) (periodicity, determinants). We also used support-vector machines (SVM) with Gaussian radial basis function (RBF) as implemented by Pedregosa et al. (2011) as baseline. RBF-SVM can be regarded as neural network with one hidden layer. Unless otherwise stated, we report the best performance across SVM parameters C' {1, 10, 100} and 7 (1/(n_fcatures iariae 1/ f/atuceurrs).
4.1 LEARNING DOUBLE SPIRALS
To study the extrapolation capability of the classification boundaries generated by the M-layer, we compare it to ReLU-DNN and RBF-SVM on a double spiral ("Swiss roll") classification task. The data consist of 2000 randomly generated points along two spirals consisting of 3 rotations with coordinates in the [-10, 10 range. Uniform random noise in the [-0.5, 0.5) range is added to each input coordinate. As we are interested in extrapolation beyond the training domain, we generate a test set containing of the continuation of those spirals, wrapping around the origin one more times.

---

Under review as a conference paper at ICLR 2021
The M-layer has representation size d = 10 and matrix size 10. The DNN has two hidden layers of size 20. We report the best result out of 5 runs. On the test set, the DNN has an accuracy of 54%, while the M-Layer achieves an accuracy of 100%. The best RBF-SVM has an accuracy of 52%; similar SVM experiments can be found in (Suykens, 2001). Figure illustrates the distinctive ability of the M-layer to extrapolate beyond the training domain.
4.2 LEARNING PERIODIC DATA
We compare an M-layer, DNNs and RBF-SVM on a dataset containing the daily minimum tem- peratures in Mclbourne recorded by the Aus- tralian Bureau of Metcorology over 10 years (Datamarket.com as cited by Kaggle.com). The data is de-meaned and smoothed using 7- day windows. The training set consists of the first 9 years, while the final year is used for test- ing. We report the best result out of 5 runs.
Temperittrre (Celsius, demeuned)
training
test set
10
RBF-SVM
DNN (505 param.)
DNN (101k param.)
M- layer (405 param.)
7.5
5
2.5
2255
The M-layer has d = 2 and n = 10, resulting in a trainable parameter count of 405. During the initialization of the M-layer, a value of 10 is
-5
-7.5
subtracted from the matrix diagonal in order to make it more likely for the initial matrix to be exponentiated to have negative eigenvalues and therefore keep outputs small. The first ReI.U DNN has two hidden layers of size 20 each, re- sulting in trainable parameter count of 505, The second ReLU DNN has linear layer followed by 10 layers of size 100, resulting in a trainable parameter count of 101301, and is trained for 5000 epochs (a minimual configuration that able fit the training data). As the aim of this experi- ment is to show the ability to learn the periodicity of the input without additional engineering, we do not consider DNNs with special activation functions such sin(x).
1981
1982
1983
1984
1985
1986
1987
1988
1989
1990
Figure 3: Predictions on periodic dataset.
Figure shows the predictions made by the trained M-layer, DNNs and RBF-SVM. The M-layer is able to pick up the periodic signal, whereas the other models attempt to fit the small-scale structure of the data, but are unable to extrapolate. Moreover, the DNN with comparable number of parameters to the M-layer does not have the capacity learn the whole training set. The RBF-SVM is able to fit the training data in high detail, but, as expected, is nonetheless unable to extrapolate beyond the training domain. During training, we observe that the M-layer first searches for an approximate solution in high error range, then rapidly decreases the loss once a suitable periodic approximation found. We also observe that the results are sensitive to parameter changes. One viable solution to this, which we have explored but defer to future work, is to add an extra parameter controlling the ratio of symmetric and antisymmetric components of the matrix to train. The eigenvalues of an antisymmetric matrix are always imaginary or zero, meaning that the matrix exponential will always expresses oscillations, whereas the symmetric component will express exponential growth or decay (see also Section 3.4).
4.3 LEARNING DETERMINANTS
To demonstrate the ability of the M-layer to learn polynomials, we train M-layers, DNNs and RBF- SVMs to predict the determinant of 3 and 5 x matrices. Learning the determinant of matrix with a small network is challenging due to the size of the search space: a 5 x 5 determinant is a polynomial with 120 monomials of degree 5 in 25 variables; the generic inhomogeneous polynomial of this degree has (25++) 142506 monomials. We do not explicitly encode any special property of the determinant. From Section 3.3, we know that it is possible to express this multivariate poly- nomial with a single strictly upper-triangular and therefore nilpotent matrix; here we demonstrate that unconstrained M-layer is equally capable of learning polynomials.
The data consist of n X n matrices with entries sampled uniformly between -1 and 1. With this m!. sampling, the expected value of the square of the determinant is z". So, we expect the square of the determinant to be 2 for 3 x 3 matrix, and for 5 matrix. This means that an estimator
7

---

Under review as a conference paper at ICLR 2021
10
Iearnt determinant
learnt determinant
5
.5
mecin squared crror
0
0
20
.5
10
-.5
.5
true determinnat
true determinant
params
3121 parums
10
M-layer
RRE-VM
DNN
(a) DNN
(b) M-layer
1111
10a
training set size
Learning the 5 x 5 determinant. Each scatter plot shows 5000 points.
Figure 5:
Figure 4: Learning the determinant of a 3 matrix with various hyperparameters.
constantly guessing 0 would have a MSE x 0.2222 and ~ 0.4938 for the two matrix sizes, respectively. which provides an error baseline for the results. The size of the training set is between 21dh and 2d examples for the 3 x 3 matrices, and 2 for the x 5 matrices. The validation set is 25% of the training set size, addition to Test sets consist of 10 matrices.
The M-layer has d = 9 and n between 6 and 12 for 3 x 3 determinants, and = 24 for 5 x 5 determinants. The DNNs has or 30 neurons, for the 3 x 3 matrices, and 5 hidden layers of size 100 for the x 5 matrices.
equally-sized hidden layers, each consisting of 5, 10, 15, 20, 25
Figure shows the results of learning the determinant 3 X matrices. The M-layer architecture able to learn from fewer examples compared to the DNN. The best M-layer model learning on 2 examples achieves a MSE of 2 10-4 with 811 parameters, while the best DNN has MSE of ss 0.003 with 3121 parameters. Figure 5 shows the results of learning the determinant of 5 x matrices. M-layer with 14977 parameters outperformed a DNN with 43101 parameters, achieving a MSE of 0.279 compared to 0.0012. We observe equivalent results for the matrix permanent.
Table 1: Comparison of classification performance on image recognition datasets. We compare the M-layer fully-connected ReLU nets, RBF-SVM, NODE (Chen al., 2018), and ANODE (Dupont et al., 2019). The ReLU net has 2 hidden layers of size (100, 50) for MNIST and (50, 50) for CIFAR-10 10 and SVHN. For NODE and ANODE, we report numbers from Dupont et (2019).
ARCHITECTRRE
CONVOLUTIONAL?
PROBLEM
ACCURACY % (MEAN + S.D.)
PARAMETERS
M-LAYER
NO
MNIST
97.99 0.12
68.855
RELU (F.C.)
NO
MNIST
97.58 = 0.10
84060
RBF-SVM
NO
MNIST
98.37
NODE
YES
MNIST
96.40 = 0.50
84 000
ANODE
YES
MNIST
98.20 0.10
84 000
M-LAYER
NO
CIFAR-10
54.17 = 0.36
148
RELU (F.C.)
NO
CIFAR- 10
51.07 = 0.62
156 710
RBF-SVM
NO
CIFAR- 10
56.86
NODE
YES
CIFAR- 10
53.70 0.20
172 000
ANODE
CIFAR -0
60.60 = 0.40
172 000
M-LAYER
NO
SVHN
81.19 0.23
148 965
RELU (F.C.)
NO
SVHN
77.83 1.23
156 710
RBF-SVM
SVHN
74.45
NODE
YES
SVHN
81.00 0.60
ANODE
YES
SVHN
83.50 = 0.50
172 000
8

---

Under review as a conference paper at ICLR 2021
4.4 LEARNING IMAGE DATASETS
We train M-layers on three image classification tasks: MNIST (Lecun et al., 1998), CIFARIO (Krizhevsky, 2009), and SVHN (Neizer al., 2011). The training set is shuffied and 10% of it is used for validation. Here, for RBF-SVM only, due to the long training times, we run the clas- sification on CIFARIO and SVHN only with the parameters that performed best on MNIST. The M-layer dimensions are always d = 35 and n 30. The M-layer can be combined with convolu- tional layers for better performance, but we limit the scope of this work to the basics of the M-layer architecture. All accuracy values are averaged over at least 30 runs.
We compare the performance of the M-layer with three penerol-purpose architectures. As the M- layer is novel architecture and no additional engineering is performed to obtain the results in addition to the regularization process described above, we only compare to other generic archi- tectures. As shown in Table the M-layer outperforms the fully-connected ReLU network and the convolutional NODE architecture while using fewer parameters. Only ANODE, an improved ver- sion of NODE, outperforms the M-layer for all datasets. RBF-SVM performs better than all other models on MNIST and worse than all other models on SVHN.
To assess the performance of computing matrix exponential, recorded training times the CIFAR- 10 task on an NVIDIA V100 GPU. The training time per epoch for the M-layer is 5.28s = 0.34 using the tf. 1 na g. expm implementation and 2.8s + 0.18 using the fast approximation described in Section 3.1 with = 3, compared to 2.41s + 0.05 for the ReLU network, where all models are trained with batch size 32. The approximation using k > 3 does not affect the accuracy.
We also compute the robustness bounds of the M-layer trained on CIFAR-IO, as described in Section 3.7. We train = 20 models with din ~ 200, ||5||2 z 3, and [|M typically value between 3 and The maximum
number of examples
600
Le variation of the vector of accumulated evidences is
400
1. This results in a typical Tox bound for robustness of R 10-3 on the whole set of correctly classified CIFAR-
200
10 test samples. In comparison, an analytical approach to robustness similar to ours (Peck et al., 2017), which uses a la-yee analysis of a traditional DNN, achieves
10
10-8 10-7 max Loo perturhation
10
10 5 10
Le bounds of = 10-9. Figure 6 shows the distribution of Loo bounds obtained for the M-layer.
Figure 6: Maximum Loo perturbation on the correctly classified CIFAR-10 10 test samples that is guaranteed not to produce a misclassification.
5 CONCLUSION
This paper introduces a novel model for supervised ma- chine learning based on single matrix exponential, where the matrix to be exponentiated depends linearly on the input. The M-layer is a powerful yet mathematically elegant architecture that has universal approximator properties and that can be used to learn and extrapolate several problems that traditional DNNs have difficulty with.
The M-layer architecture has natural ability to learn input feature crosses, multivariate polynomials and periodic functions. This allows it to extrapolate learning to domains outside the training data with no additional engineering. Its performance is similar that of non-specialized architectures, with competitive training times. Thanks to its powerful but relatively simple mathematical structure, the M-layer also allows closed-form robustness bounds.
We provide TensorFlow source code for the M-layer in the from of a Keras layer.
5.1 FUTURE WORK
In future work, we will focus on larger experiments involving advanced regularization techniques inspired by the connection between the M-layer and Lie groups We will also show how the M- layer can be combined with convolutional layers into hybrid architectures for image recognition and regression.

---

Under review as a conference paper at ICLR 2021
REFERENCES
Tian Qi Chen, Yulia Rubanova, Jesse Bettencourt, and David K Duvenaud. Neural Ordinary Dif- ferential Equations. In S Bengio. H Wallach, H Larochelle, K Grauman, N Cesa-Bianchi, and R Garnett (eds.), Advances in Neural Information Processing Systems 31, pp. Associates, Inc., 2018. NeurIPS.
Curran
Datamarket.com, as cited by Kaggle.com. www.iiimmmtteeeeeeeeeeeeeeeeeeee
com/ paulbrabbar/ htn-melbournee Accessed: 2010-10-01.
Emilien Dupont, Arnaud Doucet, and Yee Whye Teh. Augmented neural odes. In H. Wallach, H. Larocellee, A. Beygelzimer, Alché Buc, E. Fox, and R. Garnett (eds.), Advances in Neural Information Processing Systems 32, pp- 3144-3144. Curran Associates, Inc., 2019.
Melih Engin, Lei Wang, Luping Zhou, and Xinwang Liu. Deepkspd: Learning kernel-matrxx matrrix-bxeedd spd representation for fine-grained image recognition. In Proceedings of the European Conference on Computer Vision (ECCV), September 2018.
Robert Gilmore. Lie groups, physics, and geometry: an introduction for physicists, engineers and chemists. Cambridge University Press, 2008
Robert Gilmore. Lie groups, Lie algebras, and some of their applications. Courier Corporation, 2012.
Xavier Glorot and Yoshua Bengio. Understanding the difficulty of training deep feedforward neural networks. In Proceedings of the thirteenth international conference on artificial intelligence and statistics, pp. 249-256, 2010.
lan J. Goodfellow, David Warde-Farley, Mehdi Mirza, Aaron Courville, and Yoshua Bengio. Maxout networks. In 30th International Conference on Machine Learning, 2013, number PART 3, pp. 2356-2364, 2013.
Nicholas Higham. The Scaling and Squaring Method for the Matrix Exponential Revisited, SIAM Journal on Matrix Analysis and Applications, 26(4): 1779-199,, 2005. ISSN 0895-4798. doi: 10.1137/040611011
Piotr Koniusz, Hongguang Zhang, and Fatih Porikli. A deeper look power normalizations. In Proceedings of the IEFE Conference Computer Vision and Pattern Recognition (CVPR), June 2018.
Alex Krizhevsky. Learning Multiple Layers of Features from Tiny Images. Technical report, 2009.
Y. Lecun, L. Bottou, Y. Bengio, and Haffner. Gradient-based learning applied to document recognition. Proceedings of the IEEE, 86(11)2278-2324, 1998. ISSN 00189219. doi: 10.1109/5.72679.
Peihua Li, Jiangtao Xie, Qilong Wang, and Wangmeng Zuo. Is second-order information helpful for large-scale visual recognition? In Proceedings of the IEEE International Conference on Computer Vision (ICCV), Oct 2017.
Tsung-Y Yu Lin and Subhransu Maji. Improved bilinear pooling with enns. arXiv preprint arXiv: 1707.66777 2017.
Tsung-Yu Yu Lin, Aruni RoyChowdhury, and Subhransu Maji. Bilinear cnn models for fine-grained visual recognition. In Proceedings of the IEEE International Conference on Computer Vision,
December 2015.
Cleve Moler and Charles Van Loan. Nineteen Dubious Ways to Compute the Exponential of Matrix, Twenty-Five Years Later. SIAM Review, 45(1):3-49, 2003. ISSN 0036- 1445. doi: 10. 0177700036145004411
Yuval Netzer, Tao Wang. Adam Coates, Alessandro Bissacco, Bo Wu, and Andrew Y Ng. Read- ing Digits in Natural Images with Unsupervised Feature Learning. NIS Workshop on Deep Learning and Unsupervised Feature Learning 2011, 2011.
10

---

Under review as a conference paper at ICLR 2021
J. Park and I. W. Sandberg. Universal Approximation Using Basial--Bass---unctii Networks. Neu- ral Computation, 3(2):246 (22):24-2577, 1991. ISSN 0899-7667. doi: 10.11122mecc..191133222222
Jonathan Peck, Joris Roels, Bart Goossens, and Yvan Saeys. Lower bounds the robustness to adversarial perturbations. In Advances in Neural Information Processing Systems, pp. 813, 2017.
Pedregosa, G. Varoquaux, A. Gramfort, V. Michel, B. Thirion, O. Grisel, M. Blondel, P. Pretten- hofer, R. Weiss, V. Dubourg, Vanderplas, A. Passos, D. Cournapeau, M. Brucher, M. Perrot, and E. Duchesnay. Scikit-learn: Machine learning in Python. Jourmal of Machine Learning Research, 12:2825-2830, 2011.
Horn Roger and R Johnson Charles. Topics in matrix analysis.
Cambridge University Press, 1994.
D. E. Rumelhart, G. E. Hinton, and James L Mcclelland. A general framework for parallel dis- tributed processing. In Parallel distributed processing: explorations in the microstructure of cognition, chapter 2. 1986. ISBN 026268053x.
Sara Sabour, Nicholas Frosst, and Geoffrey E Hinton. Dynamic Routing Between Capsules. Ad- vances Neural Information Processing Systems 30, pp. 3856-3866, 2017.
Gagandeep Singh, Timon Gehr, Markus Püschel, and Martin Vechev. Boosting robustness certifica- tion of neural networks. 2018.
Gagandeep Singh, Timon Gehr, Markus Püschel, and Martin Vechev. An abstract domain for cer- tifying neural networks. Proceedings of the ACM on Programming Languages, 3PPPL): 2019.
Burt Speelpenning. Compiling fast partial derivatives of functions given by algorithms. PhD thesis, Illinois University, Urbana (USA), 1980.
Marshall H Stone. The generalized Weierstrass approximation theorem. Mathematics Magazine, (5):237-254, 1948.
Johan A.K. Suykens. Support vector machines: A nonlinear modelling and control perspective. European Journal of Control, 7(2):31 327, 2001. ISSN 0947-3580. doi: https://doi.or 00.3166leec.73311-327 URL http: //tp:///wwwwww// wcwencedireeeceecoooooo rectt com/ science/ article/ p51/8094735800111555
Andrew Trask, Felix Hill, Scott Reed, Jack Rae, Chris Dyer, and Phil Blunsom. Neural Arith- metic Logic Units. Advances in Neural Information Processing Systems, 2018-Decem:8035- 8044, 2018. ISSN 10495258.
A APPENDIX
A.1 UNIVERSAL. APPROXIMATION THEOREM
We show that a single M-layer model that uses sufficienlly large matrix size able to express any polynomial in the input features. This is true even when restrict the matrix be exponentiated to be nilpotent or, more specifically, strictly upper triangular. So, for classification problems, M-layer architectures are superset of multivariate polynomial classifiers, where matrix size constrains the complexity of the polynomial.
Theorem (Expressibility polynomials). Given a polynomial p(z,...., v.) in n variables, we can choose weight tensors for the M-layer such that it computes p exactly.
Proof. The tensor contraction applied to the result matrix exponentiation can form arbitrary linear combinations, and is therefore able to compute any polynomial given matrix that contains the constituent monomials up to constant factors. Thus, it suffices to prove that we can produce arbitrary monomials in the exponentiated matrix.
11

---

Under review as a conference paper at ICLR 2021
Given monomial of degree repeated) factors of the monomial on the first upper diagonal and zeros elsewhere. Let us consider powers of U. It can be shown that all elements of U equal to 0, except for the i-th upper diagonal, and that the value in the (d 1)-th upper diagonal of Ud-1, which contains only one element, is the product of the entries of the first upper diagonal of U. This precisely the monomial m we started with. By the definition of the exponential of the matrix, exp(U) then contains (d-1)!
consider the x d matrix U, that has the d 1 (possibly
which is the monomial up constant factor.
Given a polynomial p constiting of monomials, for each S L, we form matrices U2 for the corresponding monomial m. of p, as described above. Then we build the diagonal block matrix (iaag(U,..... Ut). is clear that exp(U') diag(exp(Ui)......ep)) monomials of p in exp(U).
so we can find all
To illustrate the proof, we look a monomial m abed.
0
a
0
0
0
b
U
0
0
1
Jab jabe
27 abed
0
be I
j bed ee
exp(U) =
0
0
1
e
0
0
0
1
d
0 0
0
0
The M-layers constructed here only make use of nilpotent matrices. When using this property constraint, the size of the M-layer can be effectively halved in the implementation.
The construction from Theorem can be adapted to express not only multivariate polynomial, i.e. a function to RR, but also functions to RR which restrict to polynomial in each coordinate. This, together with the Stone-Weierstee Weierstrass theorem (Stone, 1948), implies the following:
Corollary 2. For any continuous function f: [a, 6)" + R" and any >> O, there exists an M-layer model that computes a function g such that [f In-1); g(ro, ... In-i)jl for all 0 J < m.
A. OPTIMALITY OF CONSTRUCTION
While our proof constructive, we make no claim that the size the matrix used in the proof optimal and cannot be decreased. Given multivariate polynomial of degree d with monomials, the size of the matrix we construct would be ((d ((d++1 1)2. In fact, by slightly adapting the construction, we can obtain size of matrix that is tdl + 1. Given that the total number of monomials in polynomials of n variables up degree d is ("-d), it seems likely possible construct much smaller M-layers for many polynomials. Thus, one wonders what is, for a given polynomial, the minimum matrix size to represent it with an M-layer.
As an example, we look at the determinant of 3 x 3-Matrix. If the matrix is
b d f i
then the determinant is the polynomial aei - afh - bdi + bfg + cdh - ceg. From Theorem 1, we know that it is possible express this polynomial perfectly with single M-layer. However, already
12

---

Under review as a conference paper at ICLR 2021
an M-layer of size 8 is sufficient to represent the determinant of 3 x matrix: If
0
0
f
0
0
0
0
h
0
0
0
0
0
0
0
0
2d
2
0
0
0
0
0
2gg
2:
M
(3)
0
0
0
0
0
0
3c
33
0
0
0
0
0
0
3a
0
0
0
0
0
0
0
then exp(M)
1 8
-g
0
-ef odd
bfg-bdi
-eg+-dlt
flutei
oeg mes ff beg
0 1 0
2d
2
3od faf
- 3bd
2g
3eg-3ai
3bg
0 0
0
3c
3b
Sa
0.0
the sum of exp( M) 0.7 and 11,6 exactly this determinant. The permanent of a matrix
can be computed with an almost identical matrix, by removing all minus signs.
A.2 GRADIENT BACKPROPAGATION THROUGH A MATRIX EXPONENTIAL
As shown by Speclpenning (1980) his PhD thesis, which first formalized reverse mode automatic differentiation
algorithm, can backpropagate through any computational algorithm that satisfies some basic constraints. The computation needs to be executed such way that all in- termediate quantities are retained alongside zero- -nrtiiaiill accumulators. These accumulators will, on the backward pass, collect contributions of the sensitivity of the end result with respect to the associated intermediate quantity. The total effort needed for the computation of the full gradient then is a bounded small integer multiple (typically around 5, and never above 8) of the effort needed compute the function.
A basic efficient numerical algorithm that computes accurate matrix exponentials and is backpropagation--friended in this sense uses these steps:
1. Compute some norm, such as the 2-norm, of the elements of the input matrix Mo (on finite dimensional vector spaces, all norms are topologically equivalent).
2. Find an integer i such that (1/2*) My has norm smaller than some threshold such that matrices with norm not exceeding this threshold can be computed with reasonable accuracy by evaluating some matrix rational function.
3. Compute Qo == expm( M2/ by using this matrix rational function. A simple choice here would the Taylor polynomial, expm(M) x I++ I (1 (1/33)) M + .... truncated at some moderate power. (Many numerical implementations of matrix expo- nentiation instead use Padé Approximation here.)
4. Compute the sequence'Qu. Q1 QoQu. Q2 QIQ..... .k. = Qu-1Q: (i.e. square times, retaining intermediate results). Then, expm ss Qs. due to Q2 expm(Mo/ 22)222222222222
As all intermediate quantities are retained, backpropagation is straightforward. If uses Taylor ap- proximation for small-norm matrices, the only steps involved re-scalings and matrix addition and multiplications. With the Padé approximation, which uses a matrix rational function, there also an effective matrix inversion step, i.e. the need backpropagate through linear solver.
Notably, the number matrix multiplications executed by this algorithm is input data dependent, and in the context of this work, this enables the M M-layer architecture to fit oscillatory data for a large number of oscillations. If one instead were to e.g. approximate matrix exponentiation with Taylor series alone, this would easily become infeasible. The reason for this can be most
13

---

Under review as a conference paper at ICLR 2021
easily seen when trying compute Y sin x 2T 50) by Taylor-expanding the matrix exponential expm(xz) for
-
0
G= +1
For 310, the x² /n! terms in the (alternating) Taylor series become 0.01 only late 2 843, whereas ten halving-steps and Taylor expansion to only n = 5 gives us expm(zG)1 = 0.8509 (whereas sin 310 0.8519).
This data-adjjstable-dephh property is what allows the M-Layer to successfully fit and extrapolate oscillatory signals many oscillations.
14


Reviews:
[{"review_id": "iCmjTiYUTzu", "reviewer": "AnonReviewer2", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "This work proposes a novel machine learning architecture (or M-layer) that is in the form matrix exponential. This architecture can effectively model the interaction of feature components and learn multivariate polynomial functions and periodic functions. The architecture of the M-layer is well described. Its universal approximation capability is explained and proved in the Appendix. Other properties related to periodicity, connection to Lie groups, the interpretation from the perspective of dynamical systems, and the robustness bounds are discussed. Experimental study is conducted on toy datasets and three image classification datasets to demonstrate the properties of the proposed M-layer. \n\nOverall, this is a nice piece of work. The proposed M-layer is novel, and the properties shown by this kind of architecture are interesting, especially the capability in modelling periodic functions and extrapolating data. Theoretical discussion is generally clear although some places are a bit hard to follow. Experimental study supports the claims. \n\nMeanwhile, this work could be further improved at the following points:\n\n1. The experimental study does not combine the M-layer with the convolution layers, and it only tests the M-layer on relatively simple image recognition datasets. It will be interesting to know whether the M-layer can be combined with convolutional layers and trained in an end-to-end manner. This is important to check the potential of the M-layer for the applications of visual, audio, and text data analysis, for which end-to-end training has proven to be more effective. This paper could test this setting on some larger-scale benchmark such as ImageNet if possible. \n\n2. In the experimental study, it will be interesting to see the comparison between the M-layer and Gaussian RBF kernel based SVM (RBF is mentioned in Related Work) on these toy datasets and image benchmarks. RBF-SVM can be regarded as a neural network with one hidden layer. \n\n3. This work indicates that the M-layer can better consider the cross-terms of feature components. This could be related to the recent work on second (or higher)-order visual representation, in which the correlations of the channels in a convolutional feature map are extracted and used as feature representation for the subsequent fully connected and softmax layers. This paper can discuss the potential connections with this line of research ([R1]-[R4])\n\n[R1] Tsung-Yu Lin; Aruni RoyChowdhury; Subhransu Maji, Bilinear CNN Models for Fine-Grained Visual Recognition, 2015 IEEE International Conference on Computer Vision (ICCV).\n[R2] Peihua Li; Jiangtao Xie; Qilong Wang; Wangmeng Zuo, Is Second-Order Information Helpful for Large-Scale Visual Recognition? 2017 IEEE International Conference on Computer Vision (ICCV)\n[R3] Melih Engin, Lei Wang, Luping Zhou, and Xinwang Liu, DeepKSPD: Learning Kernel-Matrix-Based SPD Representation For Fine-Grained Image Recognition, European Conference on Computer Vision ECCV 2018\n[R4] Piotr Koniusz; Hongguang Zhang; Fatih Porikli, A Deeper Look at Power Normalizations, 2018 IEEE/CVF Conference on Computer Vision and Pattern Recognition\n\n4. The capability of \"extrapolate learning\" of this M-layer is very interesting. This work uses the results in Figures 2 and 3 to demonstrate it. In addition to the double spiral data and the periodic data, is this \"extrapolate learning\" capability applicable to other more general data or patterns? Please comment.  \n\n--- Thank the authors for the detailed response. After reading the response and the comments of peer reviewers, the rating is altered as follows. ", "extra_scores": {}, "extra_text": {}}, {"review_id": "bpPA07m1Fs", "reviewer": "AnonReviewer3", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "### Summary\n\nIn this work, the authors propose the M-Layer, a neural network \"layer\" that consists in affinely transforming the input to a matrix, applying its matrix exponential, and applying another affine transformation. The parameters describing the affine transformations are the learnable parameters. \nThe authors provide a \"universal approximation result\" that shows that when using a sufficiently large latent space, the matrix exponential is able to replicate arbitrary polynomial functions of the input data, as well as periodic functions. \nIn numerical experiments on synthetic data (determinant of a matrix, swiss roll), time series data, as well as image data (MNIST, CIFAR10, SVHN), a single M-Layer compares favorably to neural ODEs as well as a very simplistic RELU network, while being outperformed by augmented neural ODE's. \nThe authors furthermore provide robustness certificates for the M-layer, and point out its connections to Lie-theory and dynamical systems interpretations of deep learning.\n\n### Decision\n\nOn the positive I think that the idea of going beyond component-wise nonlinearities is interesting and the matrix exponential, with its rich mathematical structure, makes for an intriguing candidate. The authors analyse this proposal from many different view including robustness, approximation capability, computational implementation, and empirical results. \n\nOn the negative side, I find the comparison to conventional neural networks not entirely fair. In particular, the claims of improved representational power compared to traditional neural network layers seem exagerated and vague. \nThe connections drawn to Lie theory and dynamical systems seem superficial in that they mostly restate properties of the matrix exponential instead of providing interpretations of the M layer. \nFinally, the numerical experiments are comparing a single M-layer to a very shallow dense RELU network. While I appreciate the intention of the authors to keep the experimental setup simple, comparing the performance of individual layers is a very artificial situation that does not seem to give any insight if M-layers are useful as layers for deep learning. \nIf they are not, however, this casts doubt on the usefulness and relevance of M-layers in general.\n\nWhile I do think that the idea of an M-layer is interesting and in particular the relative stability bounds are promising, I feel that the message of the paper is diminished by the unclear comparison to existing deep learning methods. \nThis makes it difficult to assess the relevance of the work which is why, for now, I tend to recommend rejection. \n\n\n### Detailed comments \n \n#### Conceptual comparison to RELU:\n\n\"While highly successful in practice, this approach also has disadvantages. In a conventional DNN, any two activations only ever get combined through summation. This means that such a network requires an increasing number of parameters to approximate more complex functions even as simple as multiplication. This approach of composing simple functions does not generalize well outside the boundaries of the training data.\" \n\nI would argue that this is misleading, since the results in later layers of a deep neural network depend on the activations in a highly nonlinear way. Furthermore, I believe there exist universal approximation results even for shallow neural networks with arbitrary width, just like the universal approximation results provided for M-layers are true only in the limit of infinitely large matrix exponentials.  They also do not provide any evidence for their claim that the lack of generalization is linked to the composing of simple functions. The authors do provide some experiments to this end on swiss roll and meteorological datasets. However these examples seem somewhat taylor-made for the exponential map approximation. Furthermore, the RELU network used is not even able to fit the training data well, which suggests this issue might be more due to a lack of model capacity.\n\n#### Lie Algebras and dynamical systems\n\nMy understanding is that these sections mostly recapitulate how the matrix exponential appears in lie theory and the solution of linear ODE. However, I don't think they provide a compelling interpretation of the M-layer in terms of either of those points of view. In neural ODE for instance the input data to the layer can be thought of as a state of a dynamical system that is determined by the model weights. \nIn the M-layer, instead, the input becomes part of the system matrix, then the propagator of this system matrix (the exponential map) is computed, and finally this propagator is applied to initial conditions given by another set of weights by the M-layer. \nI cannot think of a meaningful interpretation and I don't think the authors have provided one, either.\n\n#### Numerical experiments\n\nReiterating what I wrote earlier, I really think that the performance of deep neural networks constructed with M-Layers needs to be compared to that of conventional deep neural networks to establish the relevance of the M-layer as a neural network layer. \nIn a nutshell: if the differences between M-layers and RELUs were to disappear when working with deeper networks, then why should one still use the M-layer. \n\n=====================================================================================\n\nAfter reading the rebuttal and the other reviews I still lean towards rejection, the main reason being that I am not convinced at this point that exponential nonlinearities are a direction worth pursuing.\nI furthermore find the passages on dynamical systems and lie groups to be still lacking, for the same reasons as detailed above.", "extra_scores": {}, "extra_text": {}}, {"review_id": "VoEy1ZiIxO", "reviewer": "AnonReviewer1", "summary": "", "questions": "", "limitations": "", "rating": 4, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "This paper explores matrix exponentiation as an alternative non-linearity in a neural network. The key idea is to compute an affine transform of the inputs, there by generating an nxn feature map, followed by applying a matrix exponential to this feature map, that is subsequently used for classification. The paper provides several scenarios where matrix exponential could be an interesting non-linearity to use. Experiments are provided on synthetic examples, and standard image recognition benchmarks in a limited setting and show some promise. \n\nPros:\n1. The idea of using matrix exponential as a non-linearity is an interesting idea that is perhaps not explored in the context of neural networks so far. \n2. The paper provides several contexts in which the matrix exponential can be a good non-linearity to adopt.\n3. Experiments on out-of-sample extrapolation shows promising results.\n\nCons:\n1. I think the key factor missing in the paper is perhaps a lack of a solid motivation for using the matrix exponential. What does such an operator capture intuitively in a neural network that is not possible with prior activations, for example, tanh, elu, etc.?\n\n2. Note that the matrix exponentiation operator has a complexity of O(n^1.5) for a feature matrix of nxn, and demands significant computational expense. For large n, which is what is typically used in deep networks, this operation could be impossible. It is unclear if this operation can be done on a GPU and if so to what approximation accuracy? \n\n3. There are other matrix operations that have been attempted before; for example, matrix log-maps, which can also offer non-linearities, and have been used in bilinear CNN models for fine-grained recognition tasks. The paper may include these topics in the survey and contrast the technical advantages of using the exponential instead. \n\n4. The experiments are very limited, and networks that use convolutional layers are shown to significantly outperform the proposed non-linearity. \n\n5. The paper also lacks details in analyzing the feature maps produced by the non-linearity, or provide any gradient derivations, but relies on standard autograd packages for such. It is also not clear why the paper decided to use a full image as input, and why not use convolutional layers within the pipeline. The paper title is also a bit ambiguous as it is unclear what is \"intelligent\" regarding matrix exponentiation? \n\nOverall, this paper explores an interesting non-linearity, however lacks theoretical or empirical rigour in showing that the approach is practically useful. \n\n---------------------------------\nPost-Rebuttal: Thanks to the authors for the detailed response. I think the idea has potential, however (resonating with the comments of other reviewers) I still think the paper should be significantly improved for better motivation, stronger theoretical results, and more elaborate and diverse experiments that can demonstrate the effectiveness of the method. ", "extra_scores": {}, "extra_text": {}}, {"review_id": "PngaFlimuWr", "reviewer": "AnonReviewer4", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "Title: INTELLIGENT MATRIX EXPONENTIATION\n\nSummary of the paper:\n\nThe idea is to use the matrix exponential as a layer in a neural network. This basically requires transforming to a latent variable which has a rank two tensor shape (a matrix) then applying the e.g. \"expm\" function. A robustness result is given (section 3.7) and a universal approximation result (appendix A) and some experiments on mainly toy but also image classification problems, and a univariate regression on seasonal data.\n\nPros:\n\nThe paper is written in a nice easy to read way. I enjoyed reading it because it is very easy to follow, though to be fair the contribution is rather simple so this should be a given.\n\nThe robustness result is something, but it follows from a very basic result on matrix exponentials. Also, the universal approximation result is better than nothing, but also quite straightforward and constructive and simply comes from observing the matrix exponential as an arbitrary polynomial.\n\nSection 3 is nice throughout ; simple and easy to follow yet interesting.\n\nCons:\n\nI hate to be that reviewer, but the paper should go further to motivate the use of the expm function. Figure 2 is very cool indeed, but then the real world experiments are not highly motivating. The seasonality data is neat, but this can presumably be better handled by a Gaussian process with periodic covariance, anyway (which is not compared with). To motivate using expm would require far more extensive experiments, even if this requires some heavy computational resources. \n\nSuggestions:\n\nFigure 1 - mention the need for a matrix layout (you seem to want to suggest a 3 by 3 matrix with your diagram).\n\nEquation 1 - is it normal to omit j and k from the l.h.s. in tensor notation ?\n\nRecommendation:\n\nThis paper is crying out for a more thorough experimental evaluation ; I am near the fence on this paper but I think it could be much more impactful with such and as such lean toward rejection.", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: This paper was reviewed by 4 reviewers who scored the paper below acceptance threshold even after the rebuttal. Reviewer 4 is concerned about motivation, Reviewer 2 rightly points out that there exist numerous works that use some form of spectral layers in a deep setting on challenging datasets - something lacking in this work. Reviewer 3 is concerned about limited discussion on lie groups and the overall benefit of expm(.). Reviewer 1 reverberates the same comments regarding insufficient experiments,  comparisons and limited motivation. We encourage authors to consider all pointers given by reviewers in any future re-submission.
Decision: Reject