Year: 2021

Abstract:
This paper presents a novel recommendation method called neuron-enhanced autoencoder based collaborative filtering (NE-AECF). The method uses an additional neural network to enhance the reconstruction capability of autoencoder. Different from the main neural network implemented in a layer-wise manner, the additional neural network is implemented in an element-wise manner. They are trained simultaneously to construct an enhanced autoencoder of which the activation function in the output layer is learned adaptively to approximate possibly complicated response functions in real data. We provide theoretical analysis for NE-AECF to investigate the generalization ability of autoencoder and deep learning in collaborative filtering. We prove that the element-wise neural network is able to reduce the upper bound of the prediction error for the unknown ratings, the data sparsity is not problematic but useful, and the prediction performance is closely related to the difference between the number of users and the number of items. Numerical results show that our NE-AECF has promising performance on a few benchmark datasets.

heading: 1 INTRODUCTION
text: Recommendation system aims to provide personalized recommendation based on various information such as user purchase records, social networks, user features, and item (or product) features. With the fast growth of E-commence, social media, and content provider, recommendation systems play more and more important roles in our daily life and have changed our life both explicitly and implicitly. In general, recommendation systems can be organized into three categories (Adomavicius & Tuzhilin, 2005; Zhang et al., 2019): content based method, collaborative filtering, and hybrid methods. The content based methods recommend similar items to a user or recommend one item to similar users, where the similarity is usually obtained from side information such as genre, occupation, and age. Collaborative filtering (CF) assumes that there exist potential correlations within both users and items, which can be implicitly used to predict unknown ratings. Hybrid methods are combinations of content based methods and CF methods (Adomavicius & Tuzhilin, 2005; Zhang et al., 2019; Su & Khoshgoftaar, 2009). CF is at the cores of many recommendation systems.
Early CF methods (Resnick et al., 1994) compute the similarity between users or items directly from the ratings to make prediction. This kind of method is also called memory based CF, which is easy to implement and has high interpretability. One limitation is that the similarity computed from the ratings is not informative owing to the high sparsity of the rating. Another line of CF is model based method (Ungar & Foster, 1998; Shani et al., 2002) that utilizes historical data to train a machine learning model such as Bayesian network (Breese et al., 1998; Miyahara et al., 2000) for recommendation. Model based methods are more effective than content based methods in learning complex hidden preference and handling the sparsity problem. Note that both content based and model based methods do not work when there are new items or users without ratings, which is known as the cold start problem. A popular strategy for solving the problem is to incorporate side information into CF methods (Adams et al., 2010; Welling et al., 2012; Zhang et al., 2017).
In the past decades, matrix factorization (Billsus & Pazzani, 1998; Mnih & Salakhutdinov, 2008; Koren et al., 2009) and matrix completion (Candès & Recht, 2009; Shamir & Shalev-Shwartz, 2014; Sun & Luo, 2015; Chen et al., 2016; Fan et al., 2019) have been extensively studied and used in CF. These methods usually exploit the potential low-rank structure of the incomplete rating matrix
via embedding items and users into a latent space of reduced dimension, where the observed ratings are approximated by the inner products of the user feature vectors and item feature vectors. The low-rankness is usually obtained by low-rank factorization (Koren et al., 2009), nuclear norm minimization (Candès & Recht, 2009), or Schatten-p quasi norm minimization (Fan et al., 2019). Particularly, Lee et al. (2016) proposed a local low-rank matrix approximation (LLORMA) that approximates the rating matrix as a weighted sum of a few low-rank matrices. LLORMA outperformed vanilla low-rank matrix completion methods in collaborative filtering, which indicates that the rating matrices in real applications may have more complicated structures rather than a single low-rank structure.
The success of neural networks and deep learning in computer vision and natural language processing inspired researchers to design neural networks for CF (Salakhutdinov et al., 2007; Dziugaite & Roy, 2015; Sedhain et al., 2015; Wu et al., 2016; Zheng et al., 2016; He et al., 2017; van den Berg et al., 2017; Fan & Cheng, 2018; Yi et al., 2020). For instance, Salakhutdinov et al. (2007) proposed a restricted Boltzmann machines (Hinton et al., 2006) based CF method called RBM-CF, which showed high performance in the Netflix challenge (Bennett & Lanning, 2007). Sedhain et al. (2015) proposed AutoRec, an autoencoder (Hinton & Salakhutdinov, 2006; Bengio et al., 2007) based CF method, which predicts unknown ratings by an encoder-decoder model x̂ = W2σ(W1x), where x denotes the incomplete ratings on one item or of one user and W1,W2 are weight matrices to optimize. Unlike RBM-CF, which is probabilistic and generative, AutoRec provides a discriminative approach. AutoRec outperformed LLORMA slightly on several benchmark datasets (Sedhain et al., 2015). In addition, adding depth is able to improve the performance of AutoRec (Sedhain et al., 2015). Inspried by Neural Autoregressive Distribution Estimator (NADE) (Larochelle & Murray, 2011) and RBM-CF (Salakhutdinov et al., 2007), Zheng et al. (2016) proposed a method called CF-NADE, in which parameters are shared between different ratings and achieved promising performance in several benchmarks. Muller et al. (2018) proposed a kernel based reparametrized neural network, in which the weight between two units is set to be a weighted kernel-function of the location vectors. The method works well in data visualization and recommendation systems. Interestingly, Yi et al. (2020) found that the expected value of the output layer of a neural network depends on the sparsity of the input data. They proposed a simple yet effective method called sparsity normalization to improve the performance of neural networks with sparse input data such as the highly incomplete rating matrices in CF.
It is worth mentioning that existing autoencoder based CF methods such as (Sedhain et al., 2015; Wu et al., 2016; Muller et al., 2018; Yi et al., 2020) use linear activation function in the output of the decoder, i.e., x̂ = WLhL−1, where WL denotes the weights of the output layer and hL−1 denotes the features given by the last hidden layer. Thus, these methods are under the assumption that the ratings are linear interactions between user features and item features, though the features can be nonlinear. Such an assumption may not be true or not optimal in real problems, especially when the data are bounded (e.g. images) or are collected by sensors (e.g. medical and chemical sensors) with nonlinear response functions. We suspect that the rating values given by users on items are from some nonlinear response functions because humans have complex emotion or decision curves (LeDoux, 2000; Baker, 2001). A naive method to incorporate nonlinear interaction is using nonlinear activation functions such as sigmoid function (with rescaling) in the output layer of the decoder, which however has much lower performance than using a linear activation function. That’s why existing autoencoder based CF methods use only linear activation function. Note that a pre-specified activation function for the output layer of the decoder may work on specific data but may be far away from the possible optimal choice. On the other hand, the theoretical analysis for autoencoder and deep learning based CF is very limited, while there have been many works on the theory of low-rank matrix completion based CF (Srebro & Shraibman, 2005; Candès & Recht, 2009; Shamir & Shalev-Shwartz, 2014; Fan et al., 2019).
Contribution. In this paper, we present a novel neural network CF method named NE-AECF, an enhanced autoencoder approach for recommendation system. NE-AECF is composed of two different neural networks, one is an autoencoder to reconstruct the incomplete rating matrix, while the other is an element-wise neural network to learn an activation function adaptively for the output layer of the autoencoder. We provide theoretical analysis for NE-AECF, which explains the superiority of our method. Specifically, we prove that the element-wise neural network can reduce the upper bound of the prediction error for the unknown ratings. We also prove that the data sparsity is not problematic but useful and the prediction performance is closely related to the difference between the number of
users and the number of items. Further, we demonstrate empirically our NE-AECF on benchmarks: MovieLen-100k and MovieLen-1M, achieving state-of-the-art results.
Notation. We use x (or X), x, and X to denote scalar, vector, and matrix respectively. We use ‖x‖ to denote the Euclidean norm of vector x, use ‖X‖F and ‖X‖2 to denote the Frobenius norm and spectral norm of matrix X respectively. The `21 norm of matrix is denoted by ‖X‖2,1 :=∑ i ‖xi‖, where xi denotes the i-th column ofX . The `∞ norm of matrix is denoted by ‖X‖∞ := maxij |Xij |. We use |S| to denote the cardinality of set S. The symbol ′ ′ denotes the Hadamard product between vectors or matrices. The symbol ′◦′ denotes function composition.

heading: 2 NEURON-ENHANCED AECF
text: Suppose we have an incomplete rating matrix X̃ = (x̃1, x̃2, . . . , x̃n) ∈ Rm×n, where m is the number of users and n is the number of items (without loss of generality). X̃ij ≥ 0 denotes the rating given by user i on item j and X̃ij = 0 indicates an unobserved rating. S denotes the set of observed ratings. We have X̃ij = Xij for all (i, j) ∈ S. Our goal is to predict the unobserved ratings Xij , (i, j) ∈ [m]× [n]\S, from X̃ . We want to learn a nonlinear function f : Rm 7→ Rm such that
n∑ i=1 ∥∥si (x̃i − f(x̃i))∥∥2 (1) is as small as possible, where si is a binary vector denoting whether the the corresponding element in xi is zero (unknown) or not. The motivation is predicting the missing entries of xi using its observed entries, though the missing entries of xi are filled by zeros before performing f . More formally, we consider the following problem
minimize f∈F ∥∥∥S (X̃ − f(X̃))∥∥∥2 F
(2)
where S = (s1, s2, . . . , sn), f is performed on each column of X̃ separately, and F denotes a hypothesis set of m to m functions. We have infinite choices for F . For example, F can be a set of functions in the form of neural network with some parameters W ∈ W , whereW denotes a set of matrices under some constraints. In this case, problem (2) defines a denoising autoencoder or stacked denoising autoencoders (Vincent et al., 2010), where the noises are introduced by filling the missing ratings with zeros.
Let f be an autoencoder with linear activation function in the output layer. Then (2) becomes
minimize W1,W2 ∥∥∥S (X̃ −W2σ(W1X̃))∥∥∥2 F + λ ( ‖W1‖2F + ‖W2‖2F ) , (3)
whereW1 ∈ Rd×m andW2 ∈ Rm×d are weights matrices to learn and λ is a nonnegative constant to control the strength of weight decay. We have omitted the bia terms for simplicity. σ denotes an activation function such as
ReLU σ(x) = max(x, 0) and Sigmoid σ(x) = 1/(1 + exp(−x)).
Note that (3) is exactly the basic model considered by Sedhain et al. (2015),Wu et al. (2016),Muller et al. (2018), and Yi et al. (2020). Once (3) is used, the following assumption is made implicitly.
Assumption 1. There exist two matricesA ∈ Rm×d andB ∈ Rd×n such that ‖S (X −AB) ‖F is small enough.
The assumption indicates that if d is much smaller than min(m,n),X can be well approximated by a low-rank matrix, which however may not always hold in real applications. Consider the following data generating model
X = h(A′B′), (4)
where h : R1 7→ R1 is an element-wise nonlinear function and A′ ∈ Rm×d, B′ ∈ Rd×n may be generated by some nonlinear functions. If the nonlinearity of h is high, X cannot be well approximated by a rank-d matrix. These analysis indicates that if the element-wise nonlinearity in generatingX is strong, (3) should use a large d to ensure a small enough training error.
The element-wise nonlinearity widely exists in real data. For example, in imaging science, the intensity of pixels are nonlinear responses of photoelectric element to spectrum. In chemical engineering, many sensors have nonlinear responses. In biomedical engineering, the dose-responses are often nonlinear curves. Hence, in collaborative filtering, the ratings may be nonlinear responses to some latent values, according to the studies on response curve in neuroscience and psychology (LeDoux, 2000; Baker, 2001).
Therefore, instead of (3), one may consider the following problem
minimize W1,W2 ∥∥∥S (X̃ − h(W2σ(W1X̃)))∥∥∥2 F + λ ( ‖W1‖2F + ‖W2‖2F ) , (5)
where h should be determined beforehand. A naive approach to determining h is choosing a bounded or partially bounded nonlinear function according to the range of the data. For example, if the data are image pixels within [0, 1], one may use Sigmoid function. If the data are nonnegative, one may use ReLU. However, such choices only considered the range of the data, which is just a small portion of the nonlinearity. Within the range, the true response functions are not necessarily linear (ReLU) or related to exponential (Sigmoid), and can be much more complicated.
As it is difficult to choose a suitable nonlinear function h in advance, we propose to learn h from the data adaptively, i.e.,
minimize W1,W2,h∈H ∥∥∥S (X̃ − h(W2σ(W1X̃)))∥∥∥2 F + λ ( ‖W1‖2F + ‖W2‖2F ) , (6)
where H denotes a hypothesis set of nonlinear functions from R1 to R1. We have different approaches to learning h. The first approach is combining various activation functions, i.e.,
hθ(z) = k∑ i θiσi(z), (7)
where σi(·) are different activation functions and θ = (θ1, . . . , θk)> are parameters to estimate. However, it is not clear whether (7) is able to approximate a wide range of nonlinear functions. The second approach is using polynomial functions, i.e.,
hθ(z) = k∑ i θzk. (8)
It is a k-order polynomial function and can well approximate any smooth functions provided that k is sufficiently large. Another approach is using a neural network, i.e.,
hΘ(z) = ΘLΘ(σΘ(ΘLΘ−1σΘ(· · ·σΘ(Θ1z) · · · ))), (9)
where Θ1 and ΘLΘ are vectors, Θ2, . . . ,ΘLΘ−1 are matrices, and σΘ is a fixed activation function. According to the universal approximation theorems (Pinkus, 1999; Sonoda & Murata, 2017; Lu et al., 2017), (9) is able to approximate any continuous functions provided that the network is wide enough or deep enough.
Since (9) is more flexible than (7) and (8) in function approximation, we propose to solve the following problem
minimize W,Θ ∥∥∥S (X̃ − hΘ(gW (X̃)))∥∥∥2 F + λW LW∑ l=1 ‖Wl‖2F + λΘ LΘ∑ l=1 ‖Θl‖2F , (10)
where W = {W1, . . . ,WLW }, Θ = {Θ1, . . . ,ΘLΘ}, and
gW (X̃) = WLW ( σW ( WLW−1σW (· · ·σW (W1X̃) · · · ) )) . (11)
In addition, we assume Wl ∈ Rdl×dl−1 , l ∈ [LW ], and Θl ∈ Rpl×pl−1 , l ∈ [LΘ]. Note that d0 = dLW = m and p0 = pLΘ = 1. Comparing (10) with (2), we see that we have replaced f by hW ◦ gΘ with Frobenius-norm constrained weight matrices. Model (10) is exactly our neuronenhanced autoencoder based collaborative filtering (NE-AECF) method. There are two different neural networks. The first one is an autoencoder defined by hΘ◦gW , which is to learn an contraction map from the incomplete rating matrix X̃ to itself or its observed entries more precisely. The second neural network is performed in an element-wise manner to learn an activation function h adaptively for the output layer of the autoencoder or stancked autoencoders. Figure 1 shows an example schematic of NE-AECF, where LW = LΘ = 2 and Z = gW (X̃).

heading: 3 GENERALIZATION ERROR BOUND OF NE-AECF
text: In this section, we analyze the capability of NE-AECF in predicting the unknown ratings ofX . Note that 1|S|‖S (X − X̂)‖ 2 F = 1 |S|‖S (X̃ − X̂)‖ 2 F := LS , where X̂ = hΘ(gW (X̃)). We have
LS = 1 |S| ∑
(i,j)∈S
`(Xij , X̂ij),
where `(Xij , X̂ij) = (Xij−X̂ij)2. Note that instead of the square loss, we may use other functions such as |Xij − X̂ij |. In the remainder of this paper, `(Xij , X̂ij) denotes a general loss. Let Sc , [m]× [n]\S, the generalization error of NE-AECF is quantified by
LSc = 1 |Sc| ∑
(i,j)∈Sc `(Xij , X̂ij),
which is a measurement of the prediction error of NE-AECF for the unknown ratings in X . We have the following generalization bound. Theorem 1. Suppose a set S of ratings of X ∈ Rm×n are observed uniformaly and randomly, which results in an incomplete matrix X̃ with unknown ratings replaced by some values such as zero. Suppose mn− |S| > |S| > 50. Let X̂ = hΘ(gW (X̃)), where hΘ is defined by (9) and gW is defined by (11). Suppose ‖Wl‖2 ≤ al, ‖Wl‖2,1 ≤ a′l, l ∈ [LW ], d̄ := max ( d1, . . . , dLW−1 ) < m, and ‖Θl‖2 ≤ bl, ‖Θl‖F ≤ b′l, l ∈ [LΘ]. Suppose the Lipschitz constants of σW and σΘ are ρ and % respectively. Suppose supi,j |` (Yij , Xij) | ≤ τ`, ` is η`-Lipschitz, and ‖X̂‖∞ ≤ µ. Then with
probability at least 1− δ over the random sampling S,
1 |Sc| ∑ (i,j)∈Sc ` ( Xij , X̂ij ) − 1 |S| ∑ (i,j)∈S ` ( Xij , X̂ij )
≤C1η`v1 ln |S| |S| + C2η`µ √ v2 ln v3 |S| + 11τ`mn√ |Sc||S| + 3τ` √ mn |Sc||S| ln 1 δ ,
(12)
where v1 = ρLW−1%LΘ−1‖X̃‖F √ lnm (∏LW l=1 al )(∑LW l=1 ( a′l al )2/3)3/2 (∏LΘ l=1 bl ) , v2 =
∑LΘ l=1 plpl−1, v3 = LΘρ LW−1%LΘ−1µ−1‖X̃‖F (∏LW l=1 al )(∏LΘ l=1 bl ) maxl
b′l bl , and C1, C2 are
some absolute constants.
First, let’s show that the bound is non-trivial. Since activation functions are often at most 1-Lipschitz, we let ρ = % = 1. Suppose a1 = · · · = aLW = 1 and b1 = · · · = bLΘ = 1. Since a′l/al ≤ dl−1, we
have (∑LW
l=1 ( a′la −1 l )2/3)3/2 ≤ L3/2W maxl a′la−1l ≤ L3/2W d̄. In addition, maxl b′l/bl ≤ maxl√pl. Then the bound in Theorem 1 becomes
LSc ≤LS + Õ
( η`L
3/2 W d̄‖X̃‖F |S|
) + Õ η`µ √∑LΘ
l=1 plpl−1 |S|  +
11τ`mn√ |Sc||S| + 3τ`
√ mn
|Sc||S| ln
1 δ .
(13)
Note that ‖X̃‖F ≤ √ mnmaxij |X̂ij |. If |S| > C3 max ( L 3/2 W d̄ √ mnmaxij |X̂ij |, ∑LΘ l=1 plpl−1 ) holds for some constant C3, the bound is non-trivial. Obviously, the condition holds if n is much larger than m and |S| is sufficiently large. A smaller d̄ leads to a tighter bound. More specifically, Theorem 1 provides the following results.

heading: A. The error bound of prediction for the unknown ratings can be reduced via including the additional (element-wise) neural network.
text: Given a fixed autoencoder, denote by L0S the training error without the element-wise neural network of NE-AECF. As discussed in Section 2, the additional neural network of NE-AECF aims to learn an activation function adaptively for the output layer of the decoder. Hence, the training error L0S can be reduced to LS . Denote L0S − LS = ∆S . In (13), %LΘ−1, ∏LΘ l=1 bl, and C2η`µ √ v2 ln v3 |S| are
introduced by the additional neural network. Denote v01 by the v1 without % LΘ−1 and ∏LΘ l=1 bl. Then we obtain L0Sc ≤ B0 while
LSc ≤ B0−∆S + (%LΘ−1 LΘ∏ l=1 bl − 1) C1η`v 0 1 ln |S| |S| + C2η`µ √ v2 ln v3 |S|︸ ︷︷ ︸
∆
.
Note that %LΘ−1 ∏LΘ l=1 bl can be very close to 1 and v2 is much smaller than |S| (provided that the additional network is not too large). Therefore, ∆ can be negative, which means the element-wise neural network can reduce the upper bound of the prediction error. On the other hand, if we do not use the additional neural network but still want reduce B to B0−∆S , we have to increase the depth or width of the main neural network, which will raise the value of v1 and hence increase the upper bound. These results verified the superiority of our NE-AECF over classical autoencoder based CF methods such as the AutoRec of (Sedhain et al., 2015).
B. Filling the unknown ratings with zeros reduces the upper bound of prediction error.
In Theorem 1, if the unknown ratings are replaced by zeros, ‖X̃‖F will decrease. Hence the bound become tighter. Specifically, (13) becomes
LSc ≤LS + Õ
( η`L
3/2 W d̄maxij |X̃ij |√
|S|
) + Õ η`µ √∑LΘ
l=1 plpl−1 |S|  +
11τ`mn√ |Sc||S| + 3τ`
√ mn
|Sc||S| ln
1 δ .
(14)
The bound is tight if √ |S| is much larger than d̄, the maximum size of the hidden layers.
C. Increasing n reduces the upper bound of prediction error.
Let the sampling rate |S|mn , ζ and network structures be fixed. Then (14) becomes
LSc ≤LS + Õ
( η`L
3/2 W maxij |X̃ij |√ ζmn/d̄2
) + Õ η`µ √∑LΘ
l=1 plpl−1 ζmn  +
11τ`√ (1− ζ)ζ2mn + 3τ`
√ 1
(1− ζ)ζmn ln
1 δ .
(15)
Therefore, when n increases, the bound becomes tighter. In addition, LSc ≤ LS when n → ∞. In real application, if there are more users than items, we need to construct a neural network such that the input is a vector of each user’s rating, where items correspond to features and users correspond samples. In other word, larger difference between the number of users and the number of items leads to tighter upper bound of the prediction error, because we can construct the autoencoder along the smaller size of the rating matrix.

heading: 4 CONNECTION WITH PREVIOUS WORK
text: The element-wise neural network of NE-AECF can be regarded as an activation function adaptively learned from the data. It is closely related to the previous work on adaptive activation functions such (Lin et al., 2013; Agostinelli et al., 2014; Hou et al., 2017; Goyal et al., 2019). For instance, Lin et al. (2013) proposed to use micro neural networks to improve the convolution operator in convolutional neural networks. Hou et al. (2017) showed that applying adaptive activation functions in the regression (second-to-last) layer of a neural network can significantly decrease the bias of the regression. Their adaptive activation function is in the form of piece-wise polynomials. We found that, empirically, in NE-AECF, the improvement given by polynomials (8) and piece-wise polynomials (Hou et al., 2017) are not significant, which may be caused by the unboundedness of polynomials.
As mentioned in Introduction, theoretical study for autoencoder and deep learning based collaborative filtering is very limited. In recent years, a few researchers have studied the generalization ability or sample complexity of deep neural networks (Bartlett et al., 2017; Neyshabur et al., 2018; Golowich et al., 2018) but their results do not apply to autoencoder based CF and our AE-NECF. Our proof for Theorem 1 has taken advantage of the result of (Bartlett et al., 2017). There are two major differences. First, our setting is collaborative filtering, in which the training samples are matrix elements. Hence we utilized the transductive Rademacher complexity proposed by (El-Yaniv & Pechyony, 2009). Second, computing the complexity of the element-wise neural network required a different method rather than that of (Bartlett et al., 2017).
Shamir & Shalev-Shwartz (2014) provided the following generalization bound for nuclear norm minimization based CF: LSc ≤ LS + O
( η`‖X̂‖∗( √ m+ √ n)
|S|
) + R, where X̂ denotes the recovered
matrix given by nuclear norm minimization and R stands for the remainder of their result (Theorem 6). Suppose the rank of X̂ is d̄. Then the term related to the nuclear norm can be as large as
O ( η` √ d̄‖X̂‖F ( √ m+ √ n)
ζmn
) or O ( η`µ √ d̄ ζ2m ) , where we have assumed m < n. According to (15),
the dominating term in our bound can be Õ ( η`µ √ d̄ ζ2m √ L3W ζd̄ n ) . Note that with the same d̄, the
training error of our NE-AECF is less than the training error of nuclear norm minimization because we are using neural networks. Now we conclude that when n is sufficiently large (compared to L3W ζd̄), our bound is tighter than that of (Shamir & Shalev-Shwartz, 2014).

heading: 5 NUMERICAL RESULTS

heading: 5.1 EXPERIMENTS ON MOVIELENS DATASETS
text: In this section, we evaluate the proposed method NE-AECF on two benchmark datasets of collaborative filtering: Movielens 100K and Movielens 1M (Harper & Konstan, 2015). These datasets contains 100 thousand (1 million) real world ratings for 1682 (3900) movies by 943 (6040) users, on a 5-stars scale. We randomly sample 90% of the known ratings as training set, leaving the remaining 10% as the test set. Among the training set, 5% are held out for hyperparameter tuning. Following Lee et al. (2016), test items or users without training observations are assigned with default rating of 3. The model performance is evaluated by the root mean squared error defined as
RMSE =
√∑ (i,j)∈Sc(Xij − X̂ij)2
|Sc| , (16)
where Sc denotes the set of test ratings. For all the experiments, we use Adam (Kingma & Ba, 2015) to minimize the objective function of NE-AECF.
We report the mean RMSE based on 50 random splits and compare our NE-AECF with the following methods: BisaMF (Koren et al., 2009), NNMF (Dziugaite & Roy, 2015), AutoRec (Sedhain et al., 2015), CF-NADE (Zheng et al., 2016), LLOMRA (Lee et al., 2016), GC-MC (van den Berg et al., 2017), Sparse FC (Muller et al., 2018), DMF+ (Yi et al., 2019), and AutoRec W/SN (Yi et al., 2020). In our NE-AECF, the main neural network has one hidden layer. The numbers of hidden units are 700 and 2000 for MovieLens-100k and MovieLens-1M respectively. For both datasets, the elementwise neural network has one hidden layer, of which the size is 200. The regularization parameters λW and λΘ were chosen from {0.01, 0.1, 1, 10, 50, 100, 200, 500}. As shown1 in Table 1, on MovieLens-100k, LLORMA outperformed all low-rank factorization methods, which indicates that a more sophisticated model (weighted sum of local low-rank models) fits the data better than a single low-rank model. As expected, our NE-AECF outperformed all baseline methods. Similarly, in Table 2, the RMSE of NE-AECF is less than those of other methods. It is worth mentioning that we have tried to increase the depth of the main neural network and the element-wise neural network, but the improvements in terms of RMSE are not significant.

heading: 5.2 EXPERIMENTS ON DOUBAN AND FLIXSTER
text: In this section, we evaluate the proposed NE-AECF on two more datasets Douban and Flixter, which are poplular benchmarks for recommendation systems. For these data, we use the preprocessed subsets and splits provided by Monti et al. (2017). These datasets both contain 3000 users and 3000
1The RMSEs of the compared methods shown in all tables are from the original papers, in which the experimental settings are the same or comparable.
items. Douban contains 136,891 ratings with density 0.0152 on a rating scale {1, 2, 3, 4, 5}. Flixster contains 26,173 ratings with density 0.0029 on a rating scale {0.5, 1, 1.5, ..., 5}, which is a bit different from the other three datasets. Among the training samples, 5% are used for hyperparameter tuning.
We report the mean RMSE on 5 repeated experiments and compare our NE-AECF with the following methods: PMF (Mnih & Salakhutdinov, 2008), GRALS (Rao et al., 2015), sRGCNN (Monti et al., 2017), GC-MC (van den Berg et al., 2017), Factorized EAE (Hartford et al., 2018) and GRAEM (Strahl et al., 2020). In our NE-AECF, the main neural network has one hidden layer, of which the size is 500, for both datasets. The structure of the element-wise neural network is the same as that used for the MovieLens datasets.
As shown in Table 3, our proposed NE-AECF outperforms other baseline method on both Douban and Flixster. Note that some of the compared methods include extra content like the side information of users and items into the model, while NE-AECF does not require extra content.

heading: 6 CONCLUSION
text: This paper presented a novel collaborative filtering method called NE-AECF. We analyzed the generalization ability for NE-AECF and showed that the element-wise neural network is useful in reducing the upper bound of the prediction error in CF. The theoretical analysis indicates that filling the unknown ratings with zeros can make the error bound tighter. It also provides a guideline to make a choice between item-based autoencoder and user-based autoencoder. These theoretical findings were further validated by the numerical results, in which our NE-AECF has state-of-the-art performance in terms of RMSE.
It is possible to incorporate side information or graph neural network into our NE-AECF, which can be a future work.

heading: A PROOF FOR THE MAIN THEOREM
text: First of all, we give the following lemmas.
Lemma 1. Let H = {H ∈ Rm×n : hij = ΘLΘσ (ΘLΘ−1(· · ·σ(Θ1zij) · · · )) ,∀(i, j) ∈ [m] × [n]; Θl ∈ Rpl×pl−1 , ‖Θl‖2 ≤ bl, ‖Θl‖F ≤ b′l,∀l ∈ [LΘ];Z ∈ Z, ‖Z‖F ≤ sz}, where the Lipschitz constant of σ is %. Suppose the covering number of Z with respect to ‖ · ‖F is upperbounded by κε and ε = ( 2%LΘ−1 ∏LΘ l=1 bl )−1 . Then the cover covering number of H with respect to ‖ · ‖F is bounded as
N (H, ‖ · ‖F , ) ≤ κε LΘ∏ l=1
( 3 √ 2%LΘ−1(LΘ + 1)sz ∏LΘ l=1 bl )plpl−1 .
Proof. See Section B.1.
Lemma 1 provides an upper bound of the covering number of the element-wise neural network. The following lemma shows an upper bound of the covering number of the main neural network.
Lemma 2 (Theorem 3.3 of Bartlett et al. (2017), reformulated). Let Z = WLW σ ( WLW−1 ( · · ·σ(W1X̃) · · · )) , where Wl ∈ Rdl+1dl , l ∈ [LW ], and
max(m, d1, . . . , dLW ) ≤ D. Denote the Lipschitz constant of σ by ρ. Suppose the reference matrices (M1, . . . ,MLW ) are given. Define
C = {FW(X̃) :W = (W1, . . . ,WLW ) , ‖Wl‖σ ≤ al, ‖Wl −Ml‖σ ≤ a′l,∀ ∈ [LW ]}.
The for any > 0,
lnN (C, , ‖ · ‖F ) ≤ ‖X̃‖2F ln 2D2
2
( ρ2(LW−1)
LW∏ l=1 a2l )( LW∑ l=1 ( a′l al )2/3)3 .
Now we can get an upper bound for the covering number of the entire neural network in NE-AECF. Lemma 3. The covering number ofHW,Θ = {HΘ ( FW(X̃) ) } with respect to ‖ · ‖F satisfies
lnN (H, ‖ · ‖F , ) ≤ 4%2(LΘ−1)ρ2(LW−1)‖X̃‖2F ln 2D2
2
( LW∏ l=1 a2l )( LΘ∏ l=1 b2l )( LW∑ l=1 ( a′l al )2/3)3
+ ( LΘ∑ l=1 plpl−1 ) ln 6LΘ% LΘ−1ρLW−1‖X̃‖F (∏LW l=1 al )(∏LΘ l=1 bl ) maxl b′l bl  . Proof. See Section B.2.
Now we can calculate the upper bound of the Rademacher complexity ofHW,Θ via using Lemma 3 and Dudley entropy integral bound. Lemma 4. Let v1 = 4ρ2(LW−1)%2(LΘ−1)‖X̃‖2F ln 2D2 (∏LW l=1 a 2 l )(∑LW l=1 ( a′l al )2/3)3 (∏LΘ l=1 b 2 l ) ,
v2 = ∑LΘ l=1 plpl−1, and v3 = 6LΘρ LW−1%LΘ−1‖X̃‖F (∏LW l=1 al )(∏LΘ l=1 bl ) maxl b ′ lb −1 l . Sup-
pose ‖hΘ ( gW (X̃) ) ‖∞ ≤ µ. The Rademacher complexity ofHW,Θ is bounded as
RS(HW,Θ) ≤ 4µ
S +
12 √ v1 + µ2v2 lnS
S + 12µ √ v2 lnµ−1v3√ S . (17)
Proof. See Section B.3.
The following lemma provides a sample complexity bound for transductive learning, which is consistent with the objective function and evaluation metric (RMSE) widely used in collaborative filtering. Lemma 5 (Corollary 1 of (El-Yaniv & Pechyony, 2009), reformulated). LetH be a fixed hypothesis set and suppose supi,j|X∈H |` (Yij , Xij) | ≤ τ`. Suppose a fixed set S of distinct indices is uniformly and randomly split to two subsets Strain and Stest, where2 |Stest| > |Strain| > 50. Then with probability at least 1− δ over the random split, we have
1 |Stest| ∑
(i,j)∈Stest
` (Yij , Xij) ≤ 1 |Strain| ∑
(i,j)∈Strain
` (Yij , Xij) + 4RS(` ◦ H)
+ 11τ` (|Strain|+ |Stest|)√
|Strain||Stest| + 3τ` √ (|Strain|+ |Stest|) |Strain||Stest| ln 1 δ
(18)
Then Theorem 1 can be proved as follows.
Proof. Accoding to the Rademacher contraction property, we have RS(` ◦ H) ≤ η`RS(H), where η` denotes the lipschitz constant of `. Using Lemma 5 with a slightly different notation and Lemma 4 where µ2v2 v1 provided that the element-wise neural network is small enough, we have
1 |Sc| ∑ (i,j)∈Sc ` ( Xij , X̂ij ) ≤ 1 |S| ∑ (i,j)∈S ` ( Xij , X̂ij )
+ η`
( 16µ
|S| + 48 √ v1 + µ2v2 ln |S| |S| + 48µ √ v2 lnµ−1v3√ |S| )
+ 11τ` (|S|+ |Sc|)√
|S||Sc| + 3τ` √ (|S|+ |Sc|) |S||Sc| ln 1 δ
≤ 1 |S| ∑ (i,j)∈S ` ( Xij , X̂ij )
+ C1η`v ′ 1 ln |S| |S| + C2η`µ
√ v2 ln v ′ 3
|S|
+ 11τ`mn√ |S||Sc| + 3τ`
√ mn
|S||Sc| ln
1 δ ,
(19)
where C1 and C2 are some fixed constants, and
v′1 = ρ (LW−1)%(LΘ−1)‖X‖F
√ lnD ( LW∏ l=1 al )( LW∑ l=1 ( a′l al )2/3)3/2(LΘ∏ l=1 bl ) ,
2We use these assumptions to simplify the theorem.
and v′3 = LΘµ −1γρLW−1%LΘ−1‖X‖F (∏LW l=1 al )(∏LΘ l=1 bl ) . Rename v′1 and v ′ 3 as v1 and v3 respectively, we finish the proof.

heading: B PROOF FOR LEMMAS
text: B.1 PROOF FOR LEMMA 1
Proof. Let SΘl := {Θl ∈ Rpl+1×pl : ‖Θl‖2 ≤ bl, ‖Θl‖F ≤ b′l}, ∀l ∈ [LΘ]. It is known that there exists an l-net S̄Θl obeying
N (SΘl , ‖ · ‖F , l) ≤ (
3b′l l
)plpl−1
such that ‖Θl − Θ̄l‖F ≤ l. We have
|hij − h̄ij | = ‖hij − h̄ij‖F = ∥∥ΘLΘσ (ΘLΘ−1(· · ·σ(Θ1zij) · · · ))− Θ̄LΘσ (Θ̄LΘ−1(· · ·σ(Θ̄1z̄ij) · · · ))∥∥F
= ∥∥ΘLΘσ (ΘLΘ−1(· · ·σ(Θ1zij) · · · ))− Θ̄LΘσ (ΘLΘ−1(· · ·σ(Θ1zij) · · · )) + Θ̄LΘσ (ΘLΘ−1(· · ·σ(Θ1zij) · · · ))− Θ̄LΘσ ( Θ̄LΘ−1(· · ·σ(Θ1zij) · · · ) ) + · · ·
+ Θ̄LΘσ ( Θ̄LΘ−1(· · ·σ(Θ1zij) · · · ) ) − Θ̄LΘσ ( Θ̄LΘ−1(· · ·σ(Θ̄1zij) · · · ) ) ∥∥ F
+ Θ̄LΘσ ( Θ̄LΘ−1(· · ·σ(Θ̄1zij) · · · ) ) − Θ̄LΘσ ( Θ̄LΘ−1(· · ·σ(Θ̄1z̄ij) · · · ) ) ∥∥ F
≤ ∥∥ΘLΘσ (ΘLΘ−1(· · ·σ(Θ1zij) · · · ))− Θ̄LΘσ (ΘLΘ−1(· · ·σ(Θ1zij) · · · ))∥∥F + ∥∥Θ̄LΘσ (ΘLΘ−1(· · ·σ(Θ1zij) · · · ))− Θ̄LΘσ (Θ̄LΘ−1(· · ·σ(Θ1zij) · · · )) ∥∥F + · · ·
+ ∥∥Θ̄LΘσ (Θ̄LΘ−1(· · ·σ(Θ1zij) · · · ))− Θ̄Lσ (Θ̄LΘ−1(· · ·σ(Θ̄1zij) · · · )) ∥∥F
+ ∥∥Θ̄LΘσ (Θ̄LΘ−1(· · ·σ(Θ̄1zij) · · · ))− Θ̄LΘσ (Θ̄LΘ−1(· · ·σ(Θ̄1z̄ij) · · · )) ∥∥F
(a) ≤%LΘ−1 |zij | ∥∥ΘLΘ − Θ̄LΘ∥∥F LΘ−1∏
l=1 ∥∥Θl∥∥2 + %LΘ−1 |zij |
∥∥Θ̄LΘ∥∥2∥∥ΘLΘ−1 − Θ̄LΘ−1∥∥F LΘ−2∏ l=1 ‖Θl‖2 + · · ·
+ %L−1 |zij | ( LΘ∏ l=2 ∥∥Θ̄l∥∥2 )∥∥Θ1 − Θ̄1∥∥F
+ %LΘ−1 ‖zij − z̄ij‖F LΘ∏ l=1 ‖Θ̄l‖2
≤%LΘ−1 szij LΘ ∏
l 6=LΘ
bl + szij LΘ−1 ∏
l 6=LΘ−1
bl + szij 1 LΘ∏ l 6=1 bl + · · ·+ ‖zij − z̄ij‖F LΘ∏ l=1 bl  . (20)
In (a), we used the facts ‖XY ‖F ≤ ‖X‖2‖Y ‖F and ‖σ(X)−σ(Y )‖F ≤ %‖X−Y ‖F recursively. It follows that
‖H − H̄‖F = √∑
ij
h2ij
(a) ≤%LΘ−1 √√√√√√∑ ij 2 s2zij  LΘ ∏
l 6=LΘ
bl + LΘ−1 ∏
l 6=LΘ−1 bl + · · ·+ 1 LΘ∏ l 6=1 bl 2 + ‖zij − ẑij‖2F ( LΘ∏ l=1 sl )2
= √ 2%LΘ−1 √√√√√ LΘ ∏ l 6=LΘ bl + LΘ−1 ∏ l 6=LΘ−1 bl + · · ·+ 1 LΘ∏ l 6=1 bl 2 ‖Z‖2F + ( LΘ∏ l=1 bl )2 ‖Z − Z̄‖2F
≤ √ 2%LΘ−1 √√√√√ LΘ ∏ l 6=LΘ bl + LΘ−1 ∏ l 6=LΘ−1 bl + · · ·+ 1 LΘ∏ l 6=1 bl 2 s2z + ( LΘ∏ l=1 bl )2 2z.
(21)
In (a), we used the fact (x + y)2 ≤ 2(x2 + y2). Let l = /( √
2LΘ)√ 2%LΘ−1sz ∏ k 6=l bk , ∀l ∈ [LΘ]. Let
z = / √ 2 √ 2%LΘ−1 ∏LΘ l=1 bl . We arrive at
‖H − H̄‖F ≤ . (22)
It means that H̄ is an -cover ofH . Then the covering number ofH is bounded as
N (H, ‖ · ‖F , )
≤N (Z, ‖ · ‖F , z) LΘ∏ l=1 N (SΘl , ‖ · ‖F , l)
≤κε LΘ∏ l=1
( 6%LΘ−1LΘszb ′ l ∏ k 6=l bk) )plpl−1
=κε LΘ∏ l=1
( 6%LΘ−1LΘszb ′ lb −1 l ∏LΘ k=1 bk )plpl−1
≤κε ( CΘ )∑LΘ l=1 plpl−1 ,
(23)
where CΘ = 6%LΘ−1LΘszγ ∏LΘ l=1 bl and γ = maxl b ′ lb −1 l . This finished the proof.
B.2 PROOF FOR LEMMA 3
Proof. It is easy to show that
sz ≥ ‖Z‖F = ∥∥∥WLW (σW (WLW−1σW (· · ·σW (W1X̃) · · · )))∥∥∥
F (24)
≥ ρLW−1‖X‖F LW∏ l=1 al. (25)
Combining Lemma 1 and Lemma 2, we have
lnN (H, ‖ · ‖F , )
≤ lnκε + ( LΘ∑ l=1 plpl−1 ) ln( CΘ )
≤4ρ 2(LW−1)%2(LΘ−1)‖X‖2F ln 2D2
2
( LW∏ l=1 a2l )( LW∑ l=1 ( a′l al )2/3)3(LΘ∏ l=1 b2l )
+ ( LΘ∑ l=1 plpl−1 ) ln 6LΘ% LΘ−1ρLW−1‖X‖F (∏LW l=1 al )(∏LΘ i=1 bl ) maxl b′l bl  . (26)
B.3 PROOF FOR LEMMA 4
Before proof, we give the following lemma, which is a variant of the Dudley entropy integral bound on Rademacher complexity.
Lemma 6 (Lemma A.5 of Bartlett et al. (2017)). LetF be a real-valued function class taking values in [0, 1] and assume that 0 ∈ F . Then
R(F|S) ≤ inf α>0 ( 4α√ S + 12 S ∫ √S α √ lnS(F|S , ‖ · ‖F , )d ) . (27)
Proof. For convenience, let
v1 = 4ρ 2(LW−1)%2(LΘ−1)‖X‖2F ln 2D2 ( LW∏ l=1 a2l )( LW∑ l=1 ( a′l al )2/3)3(LΘ∏ l=1 b2l ) ,
v2 = ∑LΘ l=1 plpl−1, and v3 = 6LΘρ LW−1%LΘ−1‖X‖F (∏LW l=1 al )(∏LΘ l=1 bl ) maxl b ′ lb −1 l . Then lnN (H, ‖ · ‖F , ) ≤ v1 2 + v2 ln (v3 ) from (26). Let µ = maxH∈H ‖H‖∞ and H̄ = {H̄ : µH̄ ∈ H}. It follows that lnN (H̄, ‖ · ‖F , ) ≤ µ−2v1 2 + v2 ln ( µ−1v3 ) .
According to Lemma 6, we have
RS(H̄) ≤ inf α>0 ( 4α√ S + 12 S ∫ √S α √ µ−2v1 2 + v2 ln ( µ−1v3 ) d ) (a)
≤ inf α>0 ( 4α√ S + 12 S ∫ √S α √ µ−2v1 + v2 2 + v2 lnµ−1v3 d ) (b) ≤ inf α>0 ( 4α√ S + 12 S ∫ √S α (√ µ−2v1 + v2 + √ v2 lnµ−1v3 ) d )
= inf α>0 ( 4α√ S + 12 S (√ µ−2v1 + v2 ln √ S α + √ v2 lnµ−1v3 (√ S − α ))) (c) ≤ 4 S + 12 √ µ−2v1 + v2 S lnS + 12(S − 1) √ v2 lnµ−1v3 S √ S
≤ 4 S
+ 12 √ µ−2v1 + v2 lnS
S +
12 √ v2 lnµ−1v3√
S .
(28)
In (a), we used the fact ln 1x ≤ 1 x2 . The inequality (b) holds according to
√ x+ y ≤ √ x + √ y. In
(c), we have let α = 1√ S , which though may not be the best choice. We arrive at
RS(HW,Θ) ≤ 4µ
S +
12 √ v1 + µ2v2 lnS
S + 12µ √ v2 lnµ−1v3√ S .

Review summaries:
In this paper, the authors study a traditional collaborative filtering problem with users' ratings, where the goal is to predict ratings that are unobserved. Specifically, the authors propose an enhanced autoencoder-based method called NE-AECF. The main idea of NE-AECF as shown in Eq(10) and Figure 1 is that it contains an additional module for each predicted rating by an autoencoder, i.e., h_\theta(g_w) in Eq(10) and the most right part near the output layer in Figure 1.

This paper proposes an imputation method that consists of two networks. One network is AE, and one network is element-wise non-linear transformation. The authors proved the proposed method has a better generalization bound for unknown ratings. Furthermore, the authors demonstrated the proposed method achieved promising results in a few datasets.  

This paper considers an autoencoder with one hidden layer. While the output layer uses a linear activation function in several papers in the literature, this paper proposes to use a nonlinear activation function in the output layer. Moreover, instead of using a fixed function (e.g., sigmoid), this non-linear activation function is learned in this paper, using a parametrized form. Generalization bounds are  derived, and experiments are conducted.