Title: LADDER POLYNOMIAL NEURAL NETWORKS

Abstract: Polynomial neural networks are essentially polynomial functions. These networks are shown to have nice theoretical properties by previous analysis, but they are hard to train when their polynomial orders are high. In this work, we devise a new type of activations and then create the Ladder Polynomial Neural Network (LPNN). The LPNN has a feedforward structure. It provides good control of its polynomial order because its order increases by 1 with each of its hidden layers. The new network can be treated with deep learning techniques such as batch normalization and dropout. As a result, it can be well trained with generic optimization algorithms regardless of its depth. In our empirical study, deep LPNN models achieve good performances in a series of regression and classification tasks.

1. INTRODUCTION:

A polynomial neural network has a unique underlying function, which is a polynomial. Well studied by mathematicians, polynomial functions have been shown to have many favorable theoretical properties. Polynomial neural networks also bridge the analysis of general neural architectures to the properties of polynomial functions. However, in practice, there is still a gap in performance between polynomial neural networks and popular feedforward neural networks. One purpose of this work is to narrow the gap and provide a new model for both theoretical study and practice.

One method of constructing polynomial neural networks is to use the quadratic function as activations in a feedforward network. Livni et al. (2014) call this specific type of networks as PNNs and show several of their theoretical properties, e.g. they are polynomial-time learnable, and they can approximate other feedforward neural networks. Kileel et al. (2019) study their functional space from the algebraic perspective. The theoretical analysis favors low-order polynomials in terms of learning time. However, PNNs cannot have an arbitrary order, as their orders grow exponentially with their layers. For learning PNN models, researchers often devise specialized learning algorithms (Livni et al., 2014; Soltani & Hegde, 2018; Du & Lee, 2018; Soltani & Hegde, 2019) . Most of these algorithms only work for PNNs with one or two hidden layers. A deep PNN will have a very large order, then the model becomes hard to train.

Polynomial learning models can also be constructed by defining a polynomial kernel over input features and network parameters (Blondel et al., 2016b; a; 2017) . Blondel et al. (2016b) show that PNNs with one hidden layer can be constructed by a second-order polynomial kernel. One can use an arbitrary order in the polynomial kernel, then the constructed learning model has the same order. Factorization machines (Rendle, 2010; Blondel et al., 2016a ) is another polynomial model constructed from the ANOVA kernel. Training of these models often needs special optimization methods. They don't have a layered structure, so it is not easy to apply deep learning techniques, such as batch normalization, to these models.

In this work, we propose a new type of polynomial neural networks, Ladder Polynomial Neural Networks (LPNNs). Particularly, we introduce a new activation function, product activation, and use it in a feedforward structure. In this new activation function, the hidden layer is multiplied to a linear transformation of the input feature. This activation always increases the polynomial order of the network by one, so we can specify an arbitrary order for the network. As a type of feedforward networks, an LPNN can be trained with standard training techniques. Further analysis shows that the function surface of LPNN has much better smooth guarantee than PNNs. As a function of learnable parameters, the network is block multiconvex with respect it weight matrices, providing opportunities of faster learning methods. The model space of LPNNs includes models constructed from polynomial kernels as special cases, indicating LPNNs are more flexible than these models.

The LPNN model is evaluated in a list of benchmark regression and classification tasks. The results indicate that the LPNN outperforms other polynomial networks in general. Its classification performances match those of a standard feedforward neural network. Our results also indicate the necessity of batch normalization and dropout when training deep LPNN models.

2. THE LADDER POLYNOMIAL NEURAL NETWORK:

We first define the general form of a feedforward neural network. Suppose the input to the neural network is a feature vector x ∈ R d0 , and denote h 0 = x. Suppose the network has L hidden layers, with each layer ∈ {1, . . . , L} takes the input h -1 and has the output h . Each layer is defined by

EQUATION

)

Here W is the weight matrix for layer , and σ(•) is the activation function. Note that superscripts of hidden vectors and weight matrices always denote layer indices, not exponents.

For notational simplicity, we omit intercept vectors for now and will add them back later.

We first define the product activation σ p (•) before we define the new network,

EQUATION

Here is the element-wise product. V is the parameter for the activation function. V has a shape of (d, d 0 ) if u has d entries. The activation is data-dependent. If u is a polynomial function of x, then σ p (u; V, x) is also a polynomial function of x with the polynomial order increased by 1. Note that the product activation is not a function of u because different x values may give the same u value but different responses from σ p (u; V, x).

We use product activations in a feedforward structure and get an LPNN. We use a different matrix V for the product activation in each layer . Suppose h L is the output of the neural network, the function of the LPNN is formally defined as h L = LP N N (x; θ),

EQUATION

Here we use θ to denote all network parameters, θ = W 1 , . . . , W L , V 1 , . . . , V L . The first hidden layer h 1 is a second order polynomial of the input, and each activation increase the order by 1, so the hidden layer h is an order ( + 1) polynomial. The entire network is an order L + 1 polynomial function of the input.

We further re-write the function with simple additions and multiplications.

EQUATION

Here e i is a one-hot vector with its i-th entry to be one. Note that W i and V i are both row vectors, so both W i h and V i x are scalars.

After expanding all product activations and taking out all summations, we write the entire network as

EQUATION

This equation further show that the polynomial order increases with the number of layers.

In this form, every monomial term has order L + 1. To include terms with different orders, we need to include intercept vectors. Suppose each layer has an intercept vector b , then the layer is defined by

EQUATION

We can re-write the input, hidden vectors, and weight matrices in the following form,

EQUATION

then we still have the previous form, ĥ = σ p ( Ŵ ĥ -1 ; V , x), and then all previous derivations apply. For notational simplicity, the following analysis continues to use notations without the intercept term.

3. ANALYSIS:



3.1. SMOOTHNESS OF THE FUNCTION SURFACE:

The smoothness of the function surface of a neural network characterizes the network's important properties, such as being robust to input perturbations. This subsection gives an upper bound of the Lipschitz constant of the LPNN function.

We first compute the first order derivative of h (x) with respective to the input x. According to (4), we have,

EQUATION

EQUATION

where diag(u) denotes a diagonal matrix with u as its diagonal.

Then we can give a bound of the first order derivative, which in turn gives an estimate of the Lipschitz constant of h as a function of x . The results are summarized in the following theorem.

Here, for the sake of simplicity, we use l 2 norm ( • ) for both vectors and matrices. Theorem 3.1. If h is defined recursively as in (1) and (2), then, for = 1, 2, . . . , L,

EQUATION

and

EQUATION

)

where ρ(•) takes the maximal singular value of its input matrix.

Proof. Our proof is based on two inequalities. The first one is from the definition of matrix norm, e.g. W h -1 ≤ ρ W h -1 . The second one is the Cauchy-Schwarz inequality, e.g.

(W h -1 ) (V x) ≤ (W h -1 ) • (V x) .

We will use mathematical induction to derive the results. For = 1,

h 1 = diag(V 1 x)W 1 x ≤ V 1 W 1 x 2 , ∇h 1 = diag(V 1 x)W 1 + diag(W 1 x)V 1 ≤ 2 V 1 W 1 x .

This means ( 10) and ( 11) hold for = 1. Now assume ( 10) and ( 11) hold for -1, then

h = diag(V x)W h -1 ≤ V x W h -1 ≤ k =1 V W x +1 ,

and

∇h = diag(V x)W (∇h -1 (x)) + diag(W h -1 (x))V ≤ V x W ∇h -1 (x) + W h -1 (x) V = ( + 1) k=1 V k W k x .

This bound indicates the relation between the smoothness and the norm of weight matrices. The result is similar to the feedforward neural networks with other activation functions whose gradients are always between -1 and 1 (Virmaux & Scaman, 2018) . The proof above is independent of the data and network parameters. We can further improve the bound by considering network weights with an approach similar to (Virmaux & Scaman, 2018 ), but we defer the further investigation to future work.

As a comparison, the same type of bound is much higher for a deep PNN model, as its order is exponential in its number of layers.

3.2. MULTICONVEXITY IN PARAMETERS:

We also study the network as a function of its parameters to understand its function surface for parameter optimization.

The network output h L is linear in a weight matrix W or V if we hold other weights as constants. We can see so if we examine the function form in (5): there is no multiplications between any two entries in W or V . With this property, the learning model is block multiconvex with respect to the model parameters. Theorem 3.2. Let h L = LP N N (x; θ) and y be the fitting target. Suppose loss(h L , y) is a convex loss function, then the training objective of an LPNN,

EQUATION

)

is block multiconvex in blocks {W 1 , . . . , W L , V 1 , . . . , V L }.

Multiconvexity enables some specialized optimization algorithms (Shen et al., 2017) . These algorithms have the potential to be combined with generic optimization algorithms and speed up the training of LPNN.

As a comparison, network weights in PNNs are exponentiated through hidden layers. We believe PNNs has a much rough surface than LPNNs.

3.3. TRAINING WITH BATCH NORMALIZATION AND DROPOUT:

Batch normalization (BN) (Ioffe & Szegedy, 2015) and dropout (Srivastava et al., 2014) are effective techniques for training deep neural networks. We can apply batch normalization and dropout to an LPNN without any modification given its layered structure. We put the BN layer after the activation per some practitioners' advice, though the original paper suggests putting it before the activation.

When LPNN has BN layers, the model in training is not a polynomial function, but the trained model with constant BN parameters are still polynomial functions. Here we want to integrate BN parameters into network weights so that previous derivations still apply. We consider one hidden layer and omit layer indices for notational simplicity.

Let [h 1 , . . . , h n ] be hidden layer values of n instances in a batch, then the batch-normalized hidden layer h i of instance i is computed by

EQUATION

)

Here the division / and the square are element-wise operations, and is a small positive number.

In the training procedure, the variance σ is a function of input instances in the batch, therefore, the hidden vector h i is not a linear function of the input x anymore. After training, the mean vector µ and the variance vector σ become constants, then the model is still a polynomial function.

Let's write the model with BN into the original function form. Let µ, σ, γ, and β are convergent values of the training procedure, then the equivalent LPNN without batch normalization has parameters defined as follows.

W = W diag(γ/(σ + )), b = -W diag(γµ/(σ + ) + β (14)

Here we need to build equivalence with the LPNN with intercept terms. These terms can be absorbed into weight matrices by ( 7).

In this result, we can see that BN changes the norm of weight matrices and then the Lipschitz constant of the network function. Based on the study by Santurkar et al. (2018) , BN can simply shrink the norms of weight matrices to avoid having steep slopes in the function surface.

Dropout can be directly applied to LPNN. In the training phase, using dropout is equivalent to removing some entries in summations of ( 5) and rescaling the summation. In the testing phase, dropout have no effect, and the trained model is just as the definition above.

3.4. RELATION WITH MODELS BASED ON POLYNOMIAL KERNELS:

In this subsection, we show that the polynomial networks constructed from polynomial kernel functions (Blondel et al., 2016b) are special cases of LPNN networks. Lemma 3.3. The polynomial kernel function P m (p, x) = (λ + p x) m with p, x ∈ R d can be written in the form of LP N N (x; θ) such that network weights in θ can be expressed by p.

Proof. Append 1 to the feature vector,

h 0 = [x , 1] . Set W 1 = [p , λ] , V = [p , λ]

for all = 1, . . . , m -1, and W ll = [1] for = 2, . . . , L, then LP N N (x; θ) is equivalent to the kernel by (5). Proof. By Appendix D.3 in (Blondel et al., 2016b) , the function of a factorization machine can be computed by

EQUATION

We can set the network as follows.

EQUATION

Here the identity matrix in W 1 and V 1 passes x to the hidden layer and get (x x) there. The hidden layer

h 1 = [(Ux) (Sx); x x]. Then we set W 2 = [1 , t ] with t = diag SU , then W 2 h 1 = F M (x; U, S).

From this analysis, we see that the LPNN is more flexible than polynomial models constructed from polynomial kernel functions. We believe the LPNN is also more flexible than the factorization machine. Compared with LPNNs, the factorization machine does not have monomials that contains the second or higher order exponential of a feature entry. We don't expect this difference brings much expressiveness to the factorization machine, though it creates difficulties for writing a factorization machine in the form of an LPNN. In this section, we evaluate the LPNN on several learning tasks. The LPNN is compared against feedforwrad network and three polynomial learning models. All models are summarized below.

Feedforward network (FF): FF uses ReLU functions as activations. We add l-2 norm regulraization to the model. The regularization weight is chosen from {1e-6, 1e-5, 1e-4, 5e-4}. When dropout is applied, the dropout rate is chosen from {0, 0.05, 0.1, 0.2, 0.4}.

Polynomial Neural Network (PNN): the model is the same as the FF except its activations are the quadratic function. It has the same hyperparameters as FF, and it is trained in the same way as FF.

Factorization Machine (FM): we use the implementation from the sklearn package (Niculae, Accessed in 2019) . The order of FM in this implementation can be 2 or 3. It add several ANOVA kernel functions (called factors) to increase model complexity. The model is also regularized by l-2 norm. The hyperparameters of FM include the order, the number of factors, and the weight of regularization. The number of factors is chosen from from {2, 4, 8, 16}, and the regularization weight is chose from the same range as FF. This implementation of FF does not have multiple outputs, so we have used one-vs-rest for multiclass classification problems.

Polynomial Kernel (PK): PK uses polynomial kernels. Other than that, PK is similar to FM. We can specify the order of the underlying polynomial function of PK. The hyperparameters of PK are the same as FM. The implementation is also from the sklearn package.

LPNN: the model is the same as the FF except its activations are product activations. Its hyperparameters are the same as FF, and it is trained in the same way as FF.

We test these models on five regression datasets (wine-quality, power-plant, kin8nm, bostonhousing, and concrete-strength) and six classification datasets (mnist, fashion-mnist, skin, sensIT, letter, covtype-b, and covtype). The mnist and fasion-mnist datasets come with the Keras package, the skin, sensIT, and covtype-b datasets are from the libSVM website, and all other datasets are from the UCI repository.

4.2. PRODUCT ACTIVATIONS:

We first examine the product activation function h = σ p (u; V, x) in a trained model. We set up an LPNN with three hidden layers and then train it on the mnist dataset. The training finishes after 20 epoches when the network has a validation accuracy of 0.984. Then we check inputs From each layer, activations of four hidden units are plot here in the same color. and responses of the activation functions at three different layers. We plot the response h i against the corresponding u i for each hidden unit i to generate a subplot. We randomly select 400 instances and plot each (h i , u i ) pair. We plot four hidden units at each of all three hidden layers and generate plots in Figure 1 .

In these results, we see that the product activation is not really a function because inputs with the same value may have different responses. It is also clear that the product activation is not linear. The behavior of the activation is versatile: the activations shown in (c) at layer 1 and (b) at layer 3 exhibit some linear behavior while the activations shown in (d) at layer 1 and (c) at layer 2 roughly approximate the quadratic function.

4.3. REGRESSION AND CLASSIFICATION:

We first apply the model to five regression tasks. We use the same data splits by Gal (Accessed in 2019) . Each dataset has 20 random splits. On each split, we run model selection through five-fold cross validation, re-train the model, and then test the model on the test set. The results are averaged over the 20 splits. For FF, PNN, and LPNN, we set three hidden layers and 50 hidden units in each hidden layer. We apply dropout and batch normalization to all the three models. We set the polynomial order to be 4 for the PK model to match the order as LPNN.

We set the order of FM to be 3. For each model, we select all hyperparameters described in the subsection above.

Table 1 tabulates RMSE of all algorithms on all datasets. On most datasets, no single algorithm outperforms others with statistical significance, but we put a * at the smallest error mean. If comparing the mean values only, LPNN performs a little worse than FF but better than other polynomial models. The PNN has very bad performances on two splits of the wine-quality dataset. We speculate that PNN is not stable when its polynomial order is high. PK has bad performances on wine-quality and boston-housing because the model does not fit the two tasksits performances are bad on most splits.

We then test these models on seven classification tasks. For each dataset, we set 30% as the test set, except for mnist and fasion-mnist datasets, which come with test sets. We do model selection for both architecture and hyperparameters on 20% of the training set. For neural networks, the number of hidden layers is chosen from {1, 2, 4}. We shrink the number of hidden units from the bottom to the top. The number of hidden units is computed by α (d out -d in ) + d out so that the number of hidden units in a middle layer is between the input dimension and the output dimension. The shrinking factor α is chose from {0.3, 0.5, 0.7, 0.8}. We also select the order for PK from {2, 3, 5} to match the order of LPNN. All other hyperparameters of a model are also selcted together with architectures.

The error rates of different models are reported in Table 2 . In general, the performance of LPNN is comparable to FF and better other polynomial models. Comparing to feedforward networks, LPNN has relatively better performance on classification tasks. We speculate the reason is that an LPNN only needs to decide discrete labels from its outputs in classification tasks while it needs to fit the exact value in regression tasks. LPNN may be not flexible enough for fitting continuous values compared to feedforward networks.

4.4. THE EFFECT OF BATCH NORMALIZATION AND DROPOUT:

In this subsection, we investigate the effect of batch normalization and dropout on LPNN. We use L ∈ {1, 2, 3, 5, 10}. For each depth, we try four configurations: using/not using batch normalization and using/not using dropout. We select other hyperparameters through model selection. We run the experiment on a regression task (concrete-strength) and a classification task (mnist).

The results are shown in Table 3 . From this result, we see that both batch normalization and dropout are needed to train a good LPNN model. On the mnist dataset, the LPNN without batch normalization has very bad performance when L = 5. Its performance drops sharply after a few epochs. This observation indicates that the LPNN without batch normalization is very unstable due to some bad optimization directions.

5. CONCLUSION:

In this paper, we have proposed LPNN, a new type of polynomial neural networks that can have an arbitrary polynomial order. The network is based on product activations. A product activation multiplies a linear transformation of the input to the hidden layer to achieve non-linear transformation of the input. The product activation always increases the polynomial order by 1. LPNN uses product activations in a feedforward structure. With modern training techniques, such as dropout and batch normalization, we can well train deep LPNN models. These models achieve competitive performances in our empirical evaluations. The LPNN has the potential to provide new insights to the theoretical study of polynomial learning models and. It is also a valuable learning method in practice.

:



Figures and Tables:

The learning models in the form of y = K k=1 π k (λ + p k x) m (Blondel et al., 2016b) are special cases of LPNN. Theorem 3.5. The second order factorization machines are special cases of LPNN.

Figure 1: Product activations of LPNN on the mnist dataset. The model has three hidden layers.From each layer, activations of four hidden units are plot here in the same color.

RMSE of different models on regression tasks

Table content: <table><tr><td colspan="2">methods wine-quality</td><td>power-plant</td><td>kin8nm</td><td colspan="2">boston-housing concrete-strength</td></tr><tr><td>FF</td><td colspan="3">0.60 ± 0.04* 4.05 ±0.17* 0.100 ± 0.002</td><td>2.82 ± 0.76*</td><td>5.10 ±0.49*</td></tr><tr><td>FM</td><td>0.73 ± 0.09</td><td>4.43 ± 0.15</td><td>0.155 ± 0.004</td><td>4.80 ± 1.14</td><td>8.52 ±0.59</td></tr><tr><td>PK</td><td>4.39 ± 5.50</td><td colspan="2">4.05 ± 0.15* 0.100 ± 0.005</td><td>41.9 ± 77.2</td><td>7.95 ±2.42</td></tr><tr><td>PNN</td><td>5.49 ± 16.5</td><td>5.83 ± 1.39</td><td>0.102 ± 0.007</td><td>4.59 ± 2.74</td><td>5.58 ±0.48</td></tr><tr><td>LPNN</td><td>0.82± 0.18</td><td>4.13 ± 0.16</td><td>0.099 ± 0.006*</td><td>4.05 ± 2.13</td><td>5.20 ±0.62</td></tr></table>

Error rates of different models on classification tasks

Table content: <table><tr><td colspan="3">methods mnist fashion-mnist</td><td>skin</td><td>sensIT</td><td>letter</td><td colspan="2">covtype-b covtype</td></tr><tr><td>FF</td><td>0.0185</td><td>0.108</td><td colspan="2">0.0313 0.176</td><td>0.096</td><td>0.113</td><td>0.146</td></tr><tr><td>FM</td><td>0.0573</td><td>0.167</td><td colspan="2">0.0439 0.260</td><td>0.546</td><td>0.208</td><td>0.575</td></tr><tr><td>PK</td><td>0.0506</td><td>0.168</td><td colspan="2">0.0039 0.225</td><td>0.248</td><td>0.191</td><td>0.494</td></tr><tr><td>PNN</td><td>0.0503</td><td>0.127</td><td colspan="2">0.0018 0.199</td><td>0.104</td><td>0.097</td><td>0.103</td></tr><tr><td>LPNN</td><td>0.0171</td><td>0.117</td><td colspan="3">0.0017 0.175 0.0729</td><td>0.117</td><td>0.140</td></tr><tr><td colspan="2">4 EXPERIMENT</td><td/><td/><td/><td/><td/></tr><tr><td colspan="2">4.1 EXPERIMENT SETUP</td><td/><td/><td/><td/><td/></tr></table>

Effect of batch normalization and dropout

Table content: <table><tr><td colspan="3">L BN and dropout only dropout</td><td>only BN</td><td>neither</td></tr><tr><td>1</td><td>7.77 ± 0.53</td><td colspan="2">7.78 ± 0.54 7.82 ± 0.55 7.76 ± 0.53</td></tr><tr><td>2</td><td>6.05 ± 0.50</td><td colspan="2">5.98 ± 0.57 6.26 ± 0.87 6.30 ± 0.96</td></tr><tr><td>3</td><td>5.20 ± 0.62</td><td colspan="2">5.12 ± 0.58 5.82 ± 1.13 6.89 ± 2.06</td></tr><tr><td>5</td><td>4.72 ± 0.66</td><td colspan="2">4.92 ± 0.78 5.46 ± 1.83 7.86 ± 3.10</td></tr><tr><td>10</td><td>5.11 ± 2.18</td><td colspan="2">4.71 ± 0.90 4.97 ± 0.97 7.49 ± 2.77</td></tr></table>

Effect of batch normalization and dropout on the mnist dataset

Table content: <table><tr><td colspan="5">L BN and dropout only dropout only BN neither</td></tr><tr><td>1</td><td>0.0191</td><td>0.0208</td><td>0.0242</td><td>0.0171</td></tr><tr><td>2</td><td>0.0191</td><td>0.0188</td><td>0.0202</td><td>0.0192</td></tr><tr><td>3</td><td>0.0170</td><td>0.0187</td><td>0.0229</td><td>0.0207</td></tr><tr><td>5</td><td>0.0207</td><td>0.7657</td><td>0.0271</td><td>0.8947</td></tr><tr><td>10</td><td>0.0230</td><td>0.0298</td><td>0.0207</td><td>0.0241</td></tr></table>

Peer review (iclr2020):
This paper proposes Ladder Polynomial Neural Networks (LPNNs) that use a new type of activation primitive -- a product activation -- in a feed-forward architecture. Unlike other polynomial architectures that grow in the order exponentially with network depth, the proposed approach gives explicit control over the order and smoothness of the network output and enables training with standard techniques. The proposed architecture is closely related to a decomposition of a k’th order multivariate polynomial function [T, x^{otimes k}] = lambda^	op (A x odot A x odot … odot A x) = lambda^	op (A x)^{odot k} where T is a symmetric tensor of polynomial coefficients and [cdot,cdot] denotes contraction. This is a shallow (one layer architecture) and sometimes referred as a Waring decomposition. In this paper, the authors propose a specific chain factorization of the polynomial (Eq 5 in the paper), where they write the factors recursively, that they name as a ladder polynomial neural network. h^ell = (W_ell h^{ell-1} odot V^{ell} x) The ladder architecture is very closely related to tensor trains (https://epubs.siam.org/doi/10.1137/090752286). I found it surprising and somewhat alarming that this literature is not being cited as these methods are also quite well known in deep learning. I like the smoothness analysis of section 3.1 -- the proof is quite easy to follow and direct. I would be quite surprised if this result would not be known in the literature in some other form but I don’t recall seeing it. On the other hand it seems to be inevitably very loose for a deep ladder network unless the network models the zero function. It would have been a valuable addition to the experimental section, if this bound would have been illustrated numerically on synthetic examples. In 3.2, The authors say that the objective is multiconvex -- I would argue that it is multilinear (apart from the regularization term, that is later introduces). The observation in 3.3, that batch-normalization or dropout can be used for this model is perhaps tangential to the main argument. These is investigated in the experimental section but I don’t see a clear conclusion. The section in 3.4 must include links to tensor decompositions beyond factorization machines. Overall, I think the paper has some merit and could be interesting for some readers, despite the fact that the contribution is not very original and the treatment could be improved in many ways.