Safe Collaborative Filtering
Authors: Riku Togashi
Tatsushi Oka
Naoto Ohsaka
Tetsuro Morimura
Venue: ICLR 2024 poster (2024)
Keywords: recommender systems
collaborative filtering
scalability
Primary area: general machine learning (i.e., none of the above)

Abstract: Excellent tail performance is crucial for modern machine learning tasks, such as algorithmic fairness, class imbalance, and risk-sensitive decision making, as it ensures the effective handling of challenging samples within a dataset. Tail performance is also a vital determinant of success for personalized recommender systems to reduce the risk of losing users with low satisfaction. This study introduces a "safe" collaborative filtering method that prioritizes recommendation quality for less-satisfied users rather than focusing on the average performance. Our approach minimizes the conditional value at risk (CVaR), which represents the average risk over the tails of users' loss. To overcome computational challenges for web-scale recommender systems, we develop a robust yet practical algorithm that extends the most scalable method, implicit alternating least squares (iALS). Empirical evaluation on real-world datasets demonstrates the excellent tail performance of our approach while maintaining competitive computational efficiency.

Published as a conference paper at ICLR 2024
SAFE COLLABORATIVE FILTERING
Riku T Tatsushi Oka" CyberAgent
Naoto Ohsaka"
Tetsuro Morimura"
rtogashi @acm. org', b Department of Economics, Keio University
mor mwwimura roremm/meeaeeetoorrreeeeeererere co.
ABSTRACT
Excellent tail performance erucial for modern machine learning tasks, such as algorithmic fairness, class imbalance, and risk-sensitive decision making, it ensures the effective handling of challenging samples within a dataset. Tail perfor- mance is also vital determinant of success for personalized recommender systems to reduce the risk of losing users with low satisfaction. This study introduces a " "safe" collaborative filtering method that prioritizes recommendation quality for less-satisfied users rather than focusing on the average performance. Our approach minimizes the conditional value at risk (CVaR), which represents the average risk over the tails of users' loss. To overcome computational challenges for web-scale recommender systems, we develop robust yet practical algorithm that extends the most scalable method, implicit alternating least squares (iALS). Empirical evaluation on real-world datasets demonstrates the excellent tail performance of our approach while maintaining competitive computational efficiency.
INTRODUCTION
Owing to the widespread implementation of collaborative filtering (CF) techniques (Hu et al., 2008;
Koren et al., 2009; Rendle al., 2022; Steck, 2019a), recommender systems have become ubiquitous in web applications and are increasingly impacting business profits. The quality of personalization is critical, particularly for users with low satisfaction, as they are essential for driving user growth, yet their disengagement poses a serious risk to the survival of the business. This is where the problem lies with conventional methods based on empirical risk minimization (ERM), which focus on performance in terms of averages with regard given to harmful errors for the tail users. In this paper, address this problem by utilizing a risk measure, conditional value at risk (CVaR) (Pflug, 2000; Rockafellar and Uryasev, 2000; 2002), which represents the average risk for tail users.
Thus far, significant research efforts have been devoted to the computational aspects of risk-averse optimization (Alexander et al., 2006; Xu and Zhang, 2009). We are also interested in the scalability of CVaR minimization for web-scale recommender systems, which must meet strict computational requirements. The challenge of risk-averse CF is scalability in two dimensions: users and items. In the context of personalized ranking, the focus has been on the scalability for the items to be ranked. Numerous studies (Rendle et al., 2009; Weimer et al., 2007: 2008a;b) have addressed the intractability of direct ranking optimization via the "score-and-sory approach, where, given user, a method predicts scores for all items and then sorts them according their scores. To attain further scalability for large item catalogue, practical methods adopt pointwise loss functions, enabling parallel optimization for users and items (Bayer et al., 2017: He et al., 2016; Hu et al., 2008). The key to these efficient methods is the separability of pointwise objectives, which allows block coordinate algorithms, such as alternating least squares (ALS). However, integrating CVaR minimization into this approach faces severe obstacle because of the non- sonoth and non- non-linear functions (i.e., check functions). which hinder the separability and parallel computation for items.
Our contribution in this paper is to devise a practical algorithm for risk-averse CF by eliminating the mismatch between CVaR minimization and personalized ranking. The paper is structured as follows. In Section 2, we show that applying CVaR minimization to matrix factorization (MF) is
"Coreeponniigg author

---

Published as a conference paper at ICLR 2024
computationally expensive. This is mainly due to the check function even with block multi- multi-convex, separable loss. In Section 3. we overcome this challenge by using a smoothing technique of quantile regression (Fernandes al., 2021; He et al., 2021; Tan et al.. 2022) and establish block coordinate solver that embarrassingly parallelizable for both users and items. The related work is discussed in Section 4, and our experiments are described in Section 5.
SETTING AND CHALLENGE
Following conventional studies (Weimer et al., 2007; 2008a;b), we view collaborative filtering ranking problem. In this setting, we have access to the set of implicit feedback (i. indicating that user has preferred item j. We denote the sets of users and items observed in S as U and V, respectively. For convenience, we also define Vi to be the set of items preferred by user i, i.e. V. = [j V (i,j) E 8), and U to be the set of users that have selected item j. i.e., Uj = fi E... S}. The aim of this problem is to learn a scoring function fe parameterized by 8, which produces a structured prediction fo( E TR having the same order of underlying i's preferences on V; we may denote the predicted score for (i.j) by fi(
Pairwise ranking optimization. Conventional studies (Lee et al., 2014; Park al., 2015; Rendle ct al., 2009) often formulate this problem as pairwise ranking optimization within the ERM framework, aiming to minimize population risk ming (. peai (i), V.)l, where (. represents the expectation over user with feedback V. Here, the loss function faair is expressed as follows:
faai (So(i) V.)
EE So(i,))),
MW
JEE
where 1 is the indicator function. Because this loss function is piece-wise constant and in- tractable, its convex upper bounds are used, such as margin hinge loss (Weimer et al., 2008b) and softplus loss (Rendle et al., 2009). However, such loss functions lead to non-separability for items; that is, loss cannot be decomposed into the sum of independent functions {gg jeey for items, pair (i)( V.) & Ejev 2n (f(((,,)))) Thus, this approach often relies on gradient-based optimizers with negative sampling techniques (Rendle and Freudenthaler, 2014), which however suffer from slow convergence (Chen et al., 2023; Yu et al., 2014).
Convex and separable upper bound. To address this non-separability issue, conventional methods utilize convex and separable loss functions, i.c.. pointwise loss (Hu et al., 2008; Zhou et al., 2008). We also adopt the following convex and separable upper bound of Eq. (1),
Σ [1 - fo(, Su > fo(i,)) jeV.
M
jeV
where 3o 2 1/IV| is a hyperparameter. The derivation is deferred to Appendix A.
We here implement the model of fa using matrix factorization (MF), which is a popular model owing to its scalability (Hu et al., 2008; Koren et al., 2009; Rendle et al., 2022). The model parameters of an MF comprise two blocks, i.e.. 0 = (U, V) where U E IR (2) and V R [V]dd are the embedding matrices of users and items, respectively. The prediction for user i is then defined by fo(i) Vu,, where 11 E d represents the i-th row of U, and the ERM objective can be expressed as follows:
min U.V
Σ iEU
(Vu, V.) + ((U,
U
where
1
Σ 2
1
Go 2 Vu,
(Vuu, v1)
(1 v[
(4)
M
2
JEV.
Here, V) = 1 ! denotes L2 regularization, with diagonal matrices A, x and A, R x representing user- and item dependeennd Tikhonov weights. Because a pointwise loss enables a highly scalable ALS algorithm owing to separability, it is widely used for scalable methods, such as implicit alternating least squares (iALS) (Hu et al., 2008). In fact,
2

---

Published as a conference paper at ICLR 2024
iALS is implemented in practical applications (Meng al., 2016) and is still known of the state-of-the toffthe-artt art methods (Rendle et al., 2022). The only difference between our loss and iALS loss is the normalization factor 1/[V; in the first term. For further details, Appendix B.
Conditional value at risk (CVaR). In contrast to ERM. which corresponds to the average case, focus more pessimistic scenario. Specifically, the objective of - aa)- CVaR minimization involves a conditional expectation,
Prob.
min (.,.)) ( V.) ((fe(i), V.) 2 L1ods,
0
where é1aa is the (1 - ax) qqaantile of the loss distribution, called value at risk (VaR) (Jorion, 2007). Since this is difficult to optimize, Rockafellar and Uryasev (2000; 2002) proposed the following reformulation:
B[g] 81
CVaR
Loss
min { a a (.,, (max(0, ((fo(i), V.) =)]
Figure 1: Illustration of CVaR
8,5
This objective is rather easy to optimize, as it is block multi-convex w.r.t. both and when the loss function is convex w.r.t. 0.
CVaR-ME. Now, we apply the above CVaR formulation to Eq. (3) as follows:
CA/ Σ max(O,,(Vu,, V2) - ) S(U.v))
1
min U.V.
E
(5)
However, optimizing this objective is still challenging. Although subgradient methods (Rockafel- lar and Uryasev, 2000) enable parallel optimization of non-smooth objectives without requiring separability, such first-order solvers are often impractical due to slow convergence. On the other hand, the obstacle in designing an efficient algorithm lies in the non-linearity of the ramp function (1)) = max(0, .r). destroys the objective's separability w.r.t. items, making impossible to solve the subproblem for each row of V in parallel, even with the separable upper bound in Eq, (4).
3 SAFER; SMOOTHING APPROACH FOR EFFICIENT RISK-AVERSS RECOMMENDER
3.1 CONVOLUTION-TYPE SMOOTHING FOR CVAR
To overcome the non-smoothness of CVaR, we utilize a smoothing technique for quantile regression (Fernandes et al., 2021; He et al., 2021; Man et al., 2022; Tan et al., 2022), which involves the integral convolution with smooth functions, called mollifiers, as studied by Friedrichs (1944); Schwartz (1951); Sobolev (1938) among others. We in- troduce the notation (-)()) = h -k(/hh). where k(-) is a kernel density function satisfying I k())du
1.6
PL a(1))
1.4
1.2
(Pi-o ka)(u). h 0.5
1.0
0.6
0.4
0.2
and h >> 0 denotes bandwidth, while K₁ represents its CDF; (u) = - (v)d We then define a smoothed check function Pi knn: IR R as follows:
0.0
-2.0
-1.5
-1.0
-0.5
0.0
0.5
1.0
1.5 2.0
Residual
Figure 2: Convolution-type smoothing.
(P1-a (u)) = (I-o(u)ka u)du.
(6)
where called the convolution operator. The resulting function P1-a kn attains strict convexity when kn is strictly convex?. We then obtain the objective of convolution-type smoothed CVaR::
min {E+ U.V.E
1 alul
E(er iel
(pr
kn)(e(V,, V.) {) +.(
U.V.E
We discuss the implementation of kn in Appendix D.
"hhe properties of smoothed check functions discussed Appendix C.4
3

---

Published as a conference paper at ICLR 2024
We hereafter denote the functional (U,V, as (1 a1-)). i block multi-convex w.r.t. U and V, then Iu, V) is also block multi-convex w.r.t. U. V, and & because (p1 kn, is convex and non-decreasing. Therefore, consider block coordinate algorithm, where cyclically update each block while fixing the other blocks to the current estimates. That is,
g(t+1) arg min V 11-a(1 ( v(O, E),
(8)
(U(++1, v(+1))) arg min I (U, V,(((+))). U.V
(9)
Because the subproblems lack a closed-form update and block separability, we develop our proposed method, Smoothing Approach For Efficient Risk-averse Recommender (SAFER2), which exploits the objective's smoothness to yield block multi-convex and separable optimization. In the following, we present the overall algorithm of 2AFER2, and Appendix contains the detailed implementation.
3.2 CONVOLUTION-TYEE SMOOTHED QUANTILE ESTIMATION
The smoothed (PI-a * k,) is twice differentiable, so natural approach to finding the solution {* is use the Newton-Raphson (NR) algorithm. At iteration < L in the (k: + 1)-th update, it estimates S as follows:
where df (t+1) ( V /(o,g/ gltt))
(t+1) S1+1 = n d( ((++1)
(++1)
(10)
Vi (U@, vl tl(t)))
where > represents the step size. The convolution-type smoothing was originally used with check functions for quantile regression, and the first apply smooth the CVaR objective. Interestingly, even if we apply the convolution-type smoothing to max(0,) .) in CVaR, optimizing for a given U and V equivalent smoothed (unconditional) quantile estimation; the proof is deferred to Appendix C.3.
Efficient loss computation. The naive computation of (Vuu, Vi) for every user is infeasible in large-scale settings owing to the score penalty || ||Vu, 12. which implies the materialization of the recovered matrix UVT with the cost of ((U) 1) We can reduce this cost by pre-computing and caching the Gram matrix VTV E Rdxd in O (IV2) and computing the loss of user i with ||Vu, = u! (V V)u, in (dL), thus leading to the overall cost of (((()) ()))d)..
Stochastic optimization for many users. Estimating & can be expensive because of the ((K))) computational cost for each NR step. particularly for many users. We can alleviate this extra cost by using stochastic algorithms
-hoorsanh and Mahoney, 2019; Xu et al., 2016), in which we randomly sample the user subset Ub C li and then estimate the direction by sample average approximation 1 (U, V) = + (a1(441) Exeu, (p1 * (<(V,, V.) &). This also allows us use backtracking line search (Armijo, 1966) determine an appropriate step size 7 while keeping computational costs in (14,, L) where L > 0 is the number of NR iterations.
3.3 PRIMAL-DUAL SPLITTING FOR PARALLEL COMPUTATION
The major obstacle to the scalability for items the row-wise coupling of stemming from non- linear composite (01 kn) ()))((()))))))))))) (((Vu,, V.)). Because (p1 k,,) is closed and convex, we can express it
as its biconjugate (p1 * kr)(r) = max; {a (p1 kn (2)}. This allows us to reformulate the original subproblem of U and V as the following saddle-poin optimization:
Σle [2 (L(Vu, V.) - { (p1 kn)" (2,) R(U. Q...v).
min max U.V
E aiZ
(11)
where z E RM represents the vector dual variables. Our algorithm alternatingly updates each block of variables (i.e., z. U and V) by solving its subproblem. As will discuss below, owing to convolution-type smoothing, the inner maximization for the dual variable z can be solved exactly and efficiently, and thus algorithm operates similarly to a block coordinate descent (Tseng, 2001) for the primal variables U and V. Through numerical experiments, we will show this algorithm converges in practice when given sufficiently large bandwidth h.
4

---

Published as a conference paper at ICLR 2024
Conjugate-free optimization.
In the step, solve the problem for 2, of each user
(++1) arg ( (((V()) V.) - g(++1) - (p1 * kn))
21
where the convex conjugate (pr * does not necessarily have a closed-form expression. Fortunately, with convolution- smoothing, we can sidestep this issue and solve each subproblem without explicitly computing (01 kn.". Let r:()) = ((V(*)) u (f) V.) g(++)) denote the residual for user
i U. Then, the first-order optimality condition for 2, is
2., * kb.) (x, (z( =0.
(12)
Exploiting the property of convex conjugate, the gradient (e * (zi) can be obtained by (p kn) (z;) = arg min (-x + (0r kx))
x (-V 1 * ky (24)) =0 Colap (I) =1 Kx(-z))²
= - 2 +1 +- V2 (p1 kz.) -K 24).
(13)
Recall that K3 is the CDF of kn and therefore has its inverse KT,² (i.e., the quantile function). Consequently, by using Eqs. (12) and (13), we obtain the closed-form solution of 2j
(f+1)
as:
Ve (p1 * (x) (z(++t)
=
((+1) = (-r))
(14)
Once we have g(++1), all 2 2((++) can be computed in parallel for users, and K₁, can also computed effortlessly when using standard kernels for kn, such as Gaussian, sigmoid, and logistic kernels
we-weighted alternating least squares. When given dual variables z, the optimization problem of U and V forms a re-weighted ERM, i.e., minu.v v [((W())->> Eieu 2i E( (Vu,, V.) + Q(U, V)). which is separable and block strongly multi-convex w.r.t. the rows of U and V. Consequently, this step is efficient to the same extent as the most efficient ALS solver of Hu et al. (2008). It is also worth noting that, as proposed in Hu al. (2008), we pre-compute Gram matrices V yTy V and UT diagMat(z) (z) U to reuse them for efficiently solve the separable subproblems in parallel (a.k.a the Gramian trick (Krichene et al., 2018; Rendle al., 2022)). It also possible to parallelize each pre- computation step by considering V = Ejev V,v, and U diagMat(z) = >ieu 2 u.u.
Tikhonov regularization. Previous observations have shown that Tikhonov weights A,, and A are critical in enhancing the ranking quality of MF (Rendle et al., 2022; Zhou et al., 2008). We hence develop regularization strategy for SAFER, based condition numbers, which characterize the numerical stability ridge-type problems. Let ( and x(2) denote the (i, il-eemeen of A,, and the (elj)-element of A,, respectively. We propose the weights for user i and item j as follows:
(0)
al (1 + BolV), (e1
A alu
Σ EN
+ 80a(41).
(15)
Here, À 0 represents the base weight that requires tuning. This strategy introduces only one hyperparameter and empirically improves final performance. See Appendix E for detailed derivation.
3.4 COMPUTATIONAL COMPLEXITY AND SCALABILITY.
The computational cost of the step is ((24)L), and for the U and steps is O ([4) [V))) which is identical to that of iALS (Hu et al., 2008). Therefore, the overall complexity per epoch is ( L + SId2 + ([Z() + (D))ds). The linear dependency on | and |> can be alleviated arbitrarily by increasing the parallel degree owing to separability. The cubic dependency on d can be circumvented by leveraging the conjugate gradient method (Tan et al., 2016) and subspace-bar block coordinate descent (Rendle et al., 2021). The sketch of the algorithm is shown in Algorithm 1. In Appendix D.4, we also discuss a variant of SAFER, to alleviate the d3 factor.
"The derivatives of smoothed check functions are provided in Appendix C.I. *For explicit expressions with some kemells, Remark 3.1 of He al. (2021).
5

---

Published as a conference paper at ICLR 2024
Algorithm 1: SAFER2 solver.
Algorithm 2: Subroutines for SAFER2
Input [Viieeu Output:
Function init ((A))
Vi U4 u. NVL (a/Vaa)
(U, V, heu, Gv) + Init ((iiie))
ViE Vj NO (o/Vd)La) Gy f L JEV V,v, vieu, & ((Vu,, V.)
for + 1 to T I do
E < Computeeeeet exi (E. (fi)iew)
for i + to [4] do
+ (1/1U) Exeu & return U, V, Veeiee Gv
2 - KE - 44)
u, - arg arg min, (in (Vu, V.) + A/ |ul3]
Function ComputeXi IE, (fi)e)):
Ge Eieu 2: u.u.
for I <- to L do Uniformly draw Uo Cu
for j < 1 to [V] do
d Ve _(()
Vj < arg min, ( (V,
[v11]]
Gv ΣEEN v.v,
Backtracking line search of
Vi E U, & + L(Vu., V.) return U
< < - y
return g
RELATED WORK
Risk-averse optimization. Coherent risk measures such as CVaR (also called expected short- fall (Acerbi and Tasche 2002)) are widely used instead of (incoherent) VaR (Artzner, 1997; Artzner et al., 1999; Jorion, 2007). Despite the desirable properties of CVaR (Pflug, 2000), its optimization is often found challenging owing to the non-smooth check function. Thus, in CVaR optimization and statistical learning with CVaR (Curi et al., 2020; Mhammedi et al., 2020; Soma and Yoshida, 2020), smoothing techniques are utilized, e.g., piece-wise quadratic smoothed plus (Alexander et al., 2006;
Xu and Zhang, 2009) and softplus (Soma and Yoshida, 2020). Distributionally robust optimization (DRO) is another prevalent approach to risk-averse learning (Rahimian and Mchrotra, 2019). Recent studies explore dual-free DRO algorithms to avoid maintaining large dual variables (Jin et al., 2021; Levy et al., 2020; Qi al., 2021). The dual-free algorithms introduce Bregman distance defined on dual variables (Lan and Zhou, 2018; Wang and Xiao, 2017). By contrast, we apply convolution-ge smoothing to CVaR for the efficient computation of and the separability via primal-dual splitting.
Quantile regression. In contrast to the least squares regression that models a conditional mean, quantile regression (QR) offers more flexibility model the entire conditional distribution (Koenker and Bassett Jr, 1978; Koenker and Hallock, 2001; Koenker et al., 2017). Our study is related to the computational aspects of QR that involve the piece-wise linear check function. Gu al. (2018) proposed a method based on the alternating direction method of multipliers (ADMM) (Boyd al., 2011) for smoothness. However, this approach is suited for our case because of non- separability for items in the penalty term of the augmented Lagrangian (See Appendix F). In line with Horowitz's 's smoothing (Horowitz, 1998), recent studies developed methods using convolution-type -type smoothing for large-scale inference (Fernandes et al., 2021). These studies realize tractable estimation using ADMM (Tan al., 2022), Frisch-Newton algorithm (He et al., 2021), and proximal gradient method (Man et al., 2022), whereas these QR methods cannot be used for our setting.
Robust recommender systems. To robustify recommender systems against the behavior of mali- cious minority users who have an incentive to bias data for economic reasons, previous studies (Mehta and Nejdl, 2008; Mehta et al., 2007) proposed methods using MF and Huber's M-estimator (Huber, 2011) 1), which is closely related to QR. There is a growing interest in controlling the performance distribution of recommender model, such as fairness-aware recommendation, where fairness towards users is essential (Do et al., 2021; Patro et al., 2020). Several recent studies have explored improving semi-worst-case performance for users. Singh et al. (2020) introduced a multi-objective optimization approach that balances reward and CVaR-based healthiness for online recommendation Wen et al. (2022) proposed method based on group DRO over given user groups instead of individual users. By contrast, to optimize the worst-case performance for individual users, Shivaswamy and Garcia-Gaan Garcia (2022) examined an adversarial learning approach, which trains two models. However, these methods do not focus on practical scalability and rely on gradient descent directly optimize their objectives.

---

Published as a conference paper at ICLR 2024
5 NUMERICAL EXPERIMENTS
5.1 EXPERIMENTAL SETUP
Datasets and evaluation protocol. We experiment with two MovieLens datasets (ML-/M and ML-20M) (Harper and Konstan, 2015) and Million Song Dataset (MSD) (Bertin-Mahieux et al., 2011). We strictly follow the standard evaluation methodology (Liang et al., 2018; Rendle et al., 2022; Steck, 2019b; Weimer et al., 2007) based on strong generalization. To generate implicit feedback datasets. we retain interactions with ratings larger than 4 of MovieLens datasets, while we use all interactions for MSD. We then consider 80% of users for training (i.e., (Viree). The remaining 10% of users in two holdout splits are used for validation and testing. During the evaluation phases, the 80% interactions of each user are disclosed to model as input to make predictions for the user, and the remaining 20% are used to compute ranking measures. We use Recall @ K (R@ K) as the quality measure of a ranked list and take the average over all testing users. We also evaluate the performance for semi-worst-case scenarios by considering the mean &@ K for worse-off users whose @ K is lower than the a-quantile among the testing users; note that setting a = 1.0 is equivalent to the average case. The mean R K with a = 1.0 may be referred to as R@ K for simplicity. For the validation measure, we used R @@20 for ML-IM and R@50 for ML-20M and MSD.
Models. We compare SAFER, to iALS (Rendle et al., 2022), ERM-MF in Eq. (3), and CVaR-MF in Eq. (5). To set strong baselines, we consider Mult- VAE (Liang et al., 2018) and recent iALS variant with Tikhonov regularization, which is competitive to the state-of- of-te-of-the-art methods (Rendle al., 2022). We do not compare pairwise ranking methods (e.g.. (Rendle al., 2009; Weston et al., 2011)) as they are known to be non-competitive on the three datasets (Liang et al., 2018; Rendle al.. 2022: Sedhain et al., 2016). We consider the instance of SAFER, with Gaussian kernel (u) = (2mh²² 1/2 exp((u² ) In all models, we initialize U and V with Gaussian noise with standard deviation /d where = 0.1 in all datasets (Rendle et al., 2022) and tune 3 and A. We set a = 0.3 in Eq. (5) and Eq. (7). For SAFER,, we also search the bandwidth and set the number of NR iterations as L = 5. For CV, we tune global learning rate for all variables in the batch subgradient method. For Mult-VAE, we tune the learning rate, batch size, and annealing parameter. For ML-20M and MSD, |4///11 = 0.1 for SAFER2 The dimensionality d of user/item embeddings is set to 32, 256, and three-layer architecture [V]] - +d + d - + [V]] because the encoded representation can be viewed as d-dimensional user embeddings with an auxiliary dimension of one for the bias term in the last fully-connected layer. During the validation and testing phases, each method optimizes an embedding (i.e.. u) of each user based on the user's 80% interactions while fixing trained V, and predicts the scores for all items. Notably, for each testing user, MF and SAFER, solve the objective of ERM-MF in Eq. (3) as each user's subproblem is separable and can be solved independently. Our source code is publicly available at https: / thub. com/ r i ktoor / /a fere- recommender. Further detailed descriptions and additional results are provided in Appendix G.
use the sub-sampled NR algorithm in the step with sampling ratio
512 for ML-IM, ML-20M,, and MSD, respectively; for Mult- VAE, we
5.2 RESULTS
Benchmark evaluation. The evaluation results for the three datasets are summarized in Table I. For the semi-worst-case scenarios (a = 0.3). SAFER; shows excellent quality in most cases, whereas baselines exhibit a decline in some cases. The quality of CVaR-MF deteriorates in all settings, emphasizing the advantage of our smoothing approach. The average-case performance (a 1.0) of SAFER, is also remarkable, which may be due to its robustness and generalization ability for limited testing samples. In fact, for the largest MSD with 50,000 testing users, SAFER, performs slightly worse than İALS in R @ 50 with a = 1.0. To break down the above results of ML-20M and MSD, Figure 3 compares the relative performance of each method with respect to iALS at each quantile level a. Here, because the instances of R @@ K with vary on different scales depending a, we show the relative performance of each method over iALS, obtained by dividing the method's s performance by that of iALS, in the y-axis of each figure. We see the clear improvement of SAFER, for smaller a while it maintains the average performance. It is particularly remarkable that SAFER, outperforms Mult- VAB in terms of R @5 for tail users in MI-20M.
K and
7

---

Published as a conference paper at ICLR 2024
Table 1: Ranking quality comparison on ML-IM, ML-20M,, and MSD.
MI.IM
MI-20M
MSD
R@
R@50
R@:20
R@:50
R@:50
R
100
R@
R0:00
R@050
R
100
R@:00
10
Models
o 1.0
a= 1.0
o==0.
00.3
o 0.3
o
1-0
o o===
1.0
a- 0.3
o= 00.33
2 1.0
o 1.0
o - 0.3
0.3
JALS
0,3450
0.4697
0.1166
0.2391
0.5263
0.6448
0.2085
0.5264
0.3590
0.4004
0.0963
0.1684
Malt-VAE
0,3329
0,4539
0.1106
0.2245
0,5313
0.6565
0.2091
0.3376
0.3482
0.4347
0,0854
0,1485
ERM-MF
0.3448
0.4700
0.1189
0,2315
0.5275
0.6441
0,2088
0.3251
0.3544
0.4540
0.0945
0.1644
CF
0.3318
0.4495
0.1061
0.2257
0.5031
0.6277
0.1975
0.3187
0.3234
0.4121
0.0781
0.1412
SAFER;
0.3517
0,4804
0.1279
0.2428
0.5308
0.6501
0.2152
0.3342
0.3585
0.4605
0.0983
0.1700
Dataset
6,168 users
2,811 movics
136.677 users
20,108 movies
371 ,355 users
41,140 songs
statisties
0.56 M interactions
9.34 M interactions
33.6 M interactions
ERM-MP/ALL
VAl--VEEAAL
CVaR-MF/IALS
SAFER//AALS
ML-20M
MSD
R@50
R0100
R050
R0100
1.05
Relative perfommnce
Relarive performance
Relarive perfomancce
1.0
Reletive performnce
1.00
1.05
7.05
0.00
3.00
1.00
0.00
0:35
0.75
0.50
0.50
Figure 3: Relative performance of each method over iALS for each quantile level.
1.0)
0.3)
R850
1.0)
R050
0.3)
Relative pafermance
1.05
Rolative reformacc
Reaaive perfor performncc
1.08
Pelatov peformance
1.06
1.04
1.00
1.04
1.03
1.05
1.01
1.02
1.02
3.01
1.00
1.00
1.00
1.00
0.99
0.98
JAL8
ERM MF
FALS
ERM-M
1A18
10.
ERM MF
ERM AI
Figure 4: Distributions of relative performances of SAFER, over ERM-MF and ÍALS on ML-IM.
Robustness. Considering the unexpected improvement of SAFER, for the average scenario on ML-IM, we conduct further analysis of the average performance for limited testing users. We repeatedly evaluate each model using the aforementioned protocol on 50 data splits independently generated from ML-IM with different random seeds. This evaluation procedure can be considered as a nested cross-vass volssaaliootion (Cawley and Talbot, 2010) with 50 outer folds and inner folds. The resulting relative performances of SAFER, over ERM-MF and iALS$ are shown in Figure 4, with each point indicating the ratio of testing measurements for a particular data split. Each measurement is obtained by taking an average of 10 models with different initialization weights on the same split. SAFER, generally exhibits a superior quality, even for the average cases. Furthermore, its advantage in the case with 0.3 highlights the robustness of the proposed approach.
Runtime comparison. The runtime per epoch of each m is shown Table 2. All experiments of MF-ed methods (i.e., iALS, ERM-MF, CVaR-MF, and SAFER:) are performed using our multi-threadee C++ implementation, orig- inally provided by Rendle et al. (2022), which utilizes Eigen to perform vector/matrix operations that support AVX instruc- tions. The reported numbers are the averaged runtime through
Table 2: Runtime per epoch.
ML-20M
MSD
Models TATS
Runtimc/cpach
Runtimc/epachi
3.16 sec
53.5 sec
Mult-VAE
9.06 sec
112.5 see
CVVR-MF
1.67
25.2
SAFER;
3.45 sec
57.0 sec
50 epochs measured using 86.4 GB RAM and Intel(R) Xeon(R) CPU @ 2.00GHz with 96 CPU cores. We implemented Mult-VAF using PyTorch and utilized an NVIDIA P100 GPU to speed up its training. We here consider a variant of CV with the Adam preconditioner (Kingma and Ba, 2014) as baseline. The runtime of SAFER, is competitive to iALS, which is the most efficient.
We omitted CVaR-MF and Mul as this protocol is costly due to their slow convergence.
We omitted FRM-MF as it nearly identical to iAL.S.
8

---

Published as a conference paper at ICLR 2024
IALS
VAE
CVARMF
CVaR-MF-MManE
SAFER;
0.5
MI-20M
1,0
03
0.6
066
D00000
0.5
0.2
0.4
0.0
13
10
10"
10"
10F
10P
10"
10P
10"
Traiting spects
Trainig spochs
0.3
MSD
Wall time
O- 1.0
Will dine
0.3
0.4
0.4
R8100
2.10
03:00
0.10
0.00
0.05
0.3
2.00
10
10
10'
10"
10P Wall time
10
10a
10'
10"
Training spects
Training apochs
Wall cime
Figure 5: Ranking quality vs. training epochs/wall time on ML-20M (top) and MSD (bottom).
0.06
0.18
1el
1
Ville
norm
10
Resiial norm
101
101
105
10%
00
20
25
10
20
Training epochs
Training epochs
Training epocts
0.85
0..05
3.0
0.34
0.5
0.16
8
0.30
0.14
0.80
0.62
0.61
10
20
30
40
50
10
20
40
50
0.10
10
20
50
Training epochs
Training spochs
Training spochs
Figure 6: Convergence profile of SAFER, with different bandwidth h on ML-20M.
CF method. CVaR-MF is faster than SAFER, and iALS in terms of runtime per epoch because is free from the cubic dependency d. However, as shown in Figure 5, CVaR-MF requires much more training epochs to obtain acceptable performance even with the Adam preconditioner; it has not converged yet even with 50 epochs while iALS and SAFER, converge with around 10 epochs. Mult-VAB is inefficient in terms of both runtime/epoch and training epochs even with GPU: it thus shows very slow convergence in terms of wall time (the right side of Figure 5). These results show that SAFER, enables efficient optimization in terms of both runtime and convergence speed.
Convergence profile. Lacking a theoretical convergence guarantee, we empirically investigate the training behavior of SAFER2 Figure 6 demonstrates the effect of bandwidth h on the convergence profile of SAFER, The top row of the figures shows the residual norms of each block at each step, while the bottom row illustrates the validation measures. We observe that setting small bandwidth (h = 0.01, 0.05) impedes convergence; in particular, the case with h = 0.01, which almost degenerates to the non-smooth CVaR, experiences fluctuation in residual norms and ranking quality. By contrast, a sufficiently large bandwidth (h = 0.18, 1e16) ensures stable convergence, whereas h- le16, which reduced ERM, does not achieve optimal performance in semi-worst- case scenarios. These results support that convolution- ttppe smoothing is vital for SAFER,
6 CONCLUSION
Towards the modernization of industrial recommender systems, where the engagement of tail users is vital for business growth, we presented a practical algorithm that ensures high-quality personalization for each individual user while maintaining scalability for real-world applications. Our algorithm called SAFER, overcomes non-smoothness and non-separability in CVaR minimization by using convolution-lype smoothing which is the essential ingredient to obtain its separable reformulation and attain scalability over two dimensions of users and items. Compared to the celebrated iALS, our SAFER, is scalable to the same extent yet exhibits superior semi-worst-cass performance without sacrificing average quality.

---

Published as a conference paper at ICLR 2024
REFERENCES
Sanjar M Abrarov and Brendan M Quine. Efficient algorithmic implementation the voigt/complex eπтor function based on exponential series approximation. Applied Mathematics and Computation. 21865):194411022 2011.
Carlo Acerbi and Dirk Tasche. Expected shortfall: a natural coherent alternative to value at risk. Economic notes, 31(2):37--388, 2002.
Siddharth Alexander, Thomas F Coleman, and Yuying Li. Minimizing evar and for portfolio of derivatives. Journal of Banking & Finance, 30(2):583-605, 2006.
Larry Armijo. Minimization of functions having lipschitz continuous first partial derivatives. Pacific Journal of mathematics, 16(1):1-3, 1966.
Philippe Artzner. Thinking coherently. Risk, 10:68-71, 1997.
Philippe Artzner, Freddy Delbaen, Jean-Mare Eber, and David Heath. Coherent measures of risk. Mathematical finance, 93):203-2228 1999.
Immanuel Bayer, Xiangnan He, Bhargav Kanagal, and Steffen Rendle. generic coordinate descent framework for learning from implicit feedback. In Proceedings of the 26th International Conference on World Wide Web, pages 1341-1350, 2017.
Thierry Bertin-Mahieux, Daniel PW Ellis, Brian Whitman, and Paul Lamere. The million song dataset. 2011.
Stephen Boyd and Lieven Vandenberghe. Convex optimization. Cambridge university press, 2004.
Stephen Boyd, Neal Parikh, Eric Chu, Borja Peleato, Jonathan Eckstein, et al. Distributed optimization and statistical learning via the alternating direction method of multipliers. Foundations and Trends@ in Machine learning, 3(1):1-122, 2011.
Gavin C Cawley and Nicola LC Talbot. On over-fitting in model selection and subsequent selection bias in performance evaluation. The Journal of Machine Learning Research, 11:2079-2107, 2010.
Chong Chen, Weizhi Ma, Min Zhang, Chenyang Wang, Yiqun Liu, and Shaoping Ma. Revisiting negative sampling vs. non-sampling in implicit recommendation. ACM Transactions on Information Systems, 41(1):1-25, 2023.
Victor Chernozhukov and Iván Fernández- Val. Subsampling inference on quantile regression pro- cesses. Sankhyā: The Indian Journal of Statistics, pages 253-276, 2005.
Sebastian Curi, Kfrr Y Levy, Stefanie Jegelka, and Andreas Krause. Adaptive sampling for stochastic
risk-averse learning. Advances in Neural Information Processing Systems, 33:1036-1047, 2020.
Virginie Do, Sam Corbett-Davies, Jamal Atif, and Nicolas Usunier. Two-sided fairness in rankings via lorenz dominance. Advances in Neural Information Processing Systems, 34:8596-8608, 2021.
Vassiliy A Epanechnikov. Non-paraeetric estimation of a multivariate probability density. Theory of Probability & Its Applications, 14(1):153-158 1969.
Marcelo Fernandes, Emmanuel Guerre, and Eduardo Horta. Smoothing quantile regressions. Journal of Business & Economic Statistics, 39(1):338-33 2021.
Kurt Outo Friedrichs. The identity of weak and strong extensions of differential operators. Trans, Am. Math. Soc., 132-11,, 1944. ISSN 0002-9947.
Walter Gautschi. Efficient computation of the complex error function. SIAM Journal on Numerical Analysis, 7(1):187-198, 1970.
Yuwen Gu, Jun Fan, Lingchen Kong, Shiqian Ma, and Hui Zou. Admm for high-dimensional sparse penalized quantile regression. Tecchometrics, 60(3):319-331, 2018.
Maxwell Harper and Joseph A Konstan. The movielens datasets: History and context. ACM Transactions Interactive Intelligent Systems, 2015.
10

---

Published as a conference paper at ICLR 2024
Xiangnan He, Hanwang Zhang, Min- Yen Kan, and Tat-Seng Chua. Fast matrix factorization for online recommendation with implicit feedback. In Proceedings of the 39th International ACM SIGIR conference Research and Development in Information Retrieval, pages 549-558, 2016.
Xuming He, Xiaoou Pan, Kean Ming Tan, and Wen- Wen-Xin Zhou. Smoothed quantile regression with large-scale inference. Journal of Econometrics, 2021.
Joel L Horowitz. Bootstrap methods for median regression models. Econometrica, pages 1327-1351, 1998.
Yifan Hu, Yehuda Koren, and Chris Volinsky. Collaborative filtering for implicit feedback datasets. In 2008 Eighth IEEE international conference on data mining, pages 263-272. Icee, 2008.
Peter Huber. Robust statistics. International encyclopedia of statistical science, pages 1248-1251. Springer, 2011.
Jikai Jin, Bohang Zhang, Haiyang Wang, and Liwei Wang. Non-convex distributionally robust optimization: Non-asymptotic analysis. Advances in Neural Information Processing Systems, 34: 2771-2782, 2021.
Philippe Jorion. Value al risk: the new benchmark for managing financial risk. The McGraw-Hii Companies, Inc., 2007.
Diederik P Kingma and Jimmy Ba. Adam: method for stochastic optimization. arXiv preprint arXiv:1412.6980, 2014.
Roger Koenker and Gilbert Bassett Jr. Regression quantiles. Econometrica: journal of the Economet- ric Society, pages 33-50, 1978.
Roger Koenker and Kevin Hallock. Quantile regression. Journal of economic perspectives, 15(4): 143-156, 2001.
Roger Koenker, Victor Chernozhukov, Xuming He, and Limin Peng. Handbook quantile regression. 2017.
Yehuda Koren, Robert Bell, and Chris Volinsky. Matrix factorization techniques for recommender systems. Computer, 4208):30-37, 2009.
Walid Krichene, Nicolas Mayoraz, Steffen Rendle, Li Zhang, Xinyang Yi, Lichan Hong, Ed Chi, and John Anderson. Efficient training on very large corpora via gramian estimation. 2018
Guanghui Lan and Yi Zhou, An optimal randomized incremental gradient method. Mathematical programming, 171 171:167-215, 2018.
Joonseok Lee, Samy Bengio, Seungyeon Kim, Guy Lebanon, and Yoram Singer. Local collaborative ranking. In Proceedings of the 23rd international conference on World wide web, pages 85-96 2014.
Daniel Levy, Yair Carmon, John Duchi, and Aaron Sidford. Large-scale methods for distributionally robust optimization. Advances in Neural Information Processing Systems, 33:8847-88888 2020.
Dawen Liang, Rahul G Krishnan, Matthew D Hoffman, and Tony Jebara. Variational autoencoders for collaborative filtering. In Proceedings of the 2018 world wide web conference, pages 689-698, 2018.
Rebeka Man, Xiaoou Pan, Kean Ming Tan, and Wen-Xin Zhou, A unified algorithm for penalized convolution smoothed quantile regression. arXiv preprint arXiv:2205.02433 2022.
Bhaskar Mehta and Wolfgang Nejdl. Attack resistant collaborative filtering. In Proceedings of the 31 annual international ACM SIGIR conference on Research and development in information retrieval, pages 75-82, 2008
Bhaskar Mehta, Thomas Hofmann, and Wolfgang Nejdl. Robust collaborative filtering. In Proceed- ings of the 2007 ACM conference on Recommender systems, pages 49-56, 2007.

---

Published as a conference paper at ICLR 2024
Xiangrui Meng. Joseph Bradley, Burak Yavuz, Evan Sparks, Shivaram Venkataraman, Davies Liu, Jeremy Freeman, DB Tsai, Manish Amde, Sean Owen, et al. Milib: Machine learning in apache spark. The Journal of Machine Learning Research, 17(1):1235-1241 2016.
Zakaria Mhammedi, Benjamin Guedj, and Robert Williamson. Pac-bayesian bound for the conditional value at risk. Advances in Neural Information Processing Systems, 33:17919-17930, 2020.
RV Mises and Hilda Dollaczek-Geiringe Praktische verfahren der gleichungsauflssngg ZAMM- Journal of Applied Mathematics and Mechanics/Zeitschrit für Angewandte Mathematik und Mechanik, 9(1)58--77, 1929.
Naoto Ohsaka and Riku Togashi. Curse "low" dimensionality recommender systems. In Proceedings of the 46th International ACM SIGIR Conference on Research and Development in Information Retrieval, page 537-547, 2023.
Dohyung Park, Joe Neeman, Jin Zhang. Sujay Sanghavi, and Inderjit Dhillon. Preference completion: Large-scale collaborative ranking from pairwise comparisons. In International Conference on Machine Learning, pages 1907-1916. PMLR, 2015.
Gourab K Patro, Arpita Biswas, Niloy Ganguly, Krishna P Gummadi, and Abhijnan Chakraborty. Fair- rec: Two-sided fairness for personalized recommendations in two-sided platforms. Proceedings of the web conference 2020, pages 1194-1204, 2020.
Georg Ch Pflug. Some remarks the value-at-risk and the conditional value-at-risk. Probabilistic constrained optimization: Methodology and applications, pages 272-281, 2000.
István Pilászy, Dávid Zibriczky, and Domonkos Tikk. Fast als-based matrix factorization for explicit and implicit feedback datasets. In Proceedings of the fourth ACM conference on Recommender systems, pages 71 71-88,7 2010.
Gert PM Poppe and Christianus MJ Wijers. More efficient computation of the complex error function. ACM Transactions on Mathematical Software (TOMS), 16(1):38 46, 1990.
Qi Qi, Zhishuai Guo, Yi Xu, Rong Jin, and Tianbao Yang. An online method for a class of distributionally robust optimization with non-convex objectives. Advances in Neural Information Processing Systems, 34: 10067-1088, 2021.
Hamed Rahimian and Sanjay Mehrotra. Distributionally robust optimization: A review. arXiv preprint wrXiv:/908.056533 2019.
Steffen Rendle and Christoph Freudenthaler. Improving pairwise learning for item recommendation from implicit feedback. In Proceedings of the 7th ACM international conference on Web search and data mining, pages 273-282, 2014.
Steffen Rendle, Christoph Freudenthaler; Zeno Gantner, and Lars Schnidt-Thieme, Bpr: Bayesian personalized ranking from implicit feedback. In Proceedings of the Twenty-Fifth Conference on Uncertainty in Artificial Intelligence, pages 452-461, 2009.
Steffen Rendle, Walid Krichene, Li Zhang, and Yehuda Koren. lals++: Speeding up matrix factoriza- tion with subspace optimization. arXiv preprint arXi::2110.14444 2021.
Steffen Rendle, Walid Krichene, Li Zhang, and Yehuda Koren. Revisiting the performance of ials on item recommendation benchmarks. In Proceedings of the 16th ACM Conference on Recommender Systems, pages
2022.
Tyrrell Rockafellar and Stanislav Uryasev. Optimization of conditional value-at-risk. Journal of risk, 2:21-42, 42, 2000.
R Tyrrell Rockafellar and Stanislav Uryasev. Conditional value-at-risk for general loss distributions. Journal of banking & finance, (447):44331444 2002.
Farbod Roosta-Khin Khorasani and Michael W Mahoney. Sub-sampled newton methods. Mathematical Programming, 174:293-326, 2019.
12

---

Published as a conference paper at ICLR 2024
Laurent Schwartz. Théorie des distributions. Hermann & Cie, 1951.
Suvash Sedhain, Aditya Menon, Scott Sanner, and Darius Braziunas. On the effectiveness of linear models for one-class collaborative filtering. In Proceedings of the AAAI Conference on Artificial Intelligence, volume 30, 2016.
Pannaga Shivaswamy and Dario Garcia-Garcia. Adversary or friend? an adversarial approach to improving recommender systems. In Proceedings of the 16th ACM Conference on Recommender Systems, pages 369-377, 2022.
Ashudeep Singh, Yoni Halpern, Nithum Thain, Konstantina Christakopoulou, E Chi, Jilin Chen, and Alex Beutel. Building healthy recommendation sequences for everyone: A safe reinforcement learning approach. In FAceTRec Workshop, 2020.
Sergei L. Sobolev. Sur un théorème d' analyse fonctionnelle. Recueil Mathématique (Matematicheski Shormik), 4(46))471 497, 1938.
Tasuku Soma and Yuichi Yoshida. Statistical learning with conditional value at risk. arXiv preprint arXiv:2002.05586, 2020.
Harald Steck. Embarrassingly shallow autoencoders for sparse data. In The World Wide Web Conference, pages 3251-3257, 2019a.
Harald Steck. Markov random fields for collaborative filtering. Advances in Neural Information Processing Systems, 32, 2019b.
Harald Steck, Maria Dimakopoulou, Nickolai Riabov, and Tony Jebara. Admm slim: Sparse recommendations for many users. In Proceedings of the 13th international conference web search and data mining, pages 555-563, 2020.
Kean Ming Tan, Lan Wang, Wen-Xin Zhou, et al. High-dimensions quantile regression: Convolution smoothing and concave regularization, Journal of the Royal Statistical Society Series 84(1): 205-233 2022.
Wei Tan, Liangliang Cao. and Liana Fong. Faster and cheaper: Parallelizing large-scale matrix factor- ization on gpus. Proceedings of the 25th ACM International Symposium on Perforrerf Parallel and Distributed Computing, pages 219-230, 2016.
Riku Togashi and Kenshi Abe. Exposure control in large-scale recommender systems. arXiv preprint 27X4v:2209.04394, 2022.
Paul Tseng. Convergence of block coordinate descent method for nondifferentiable minimization. Journal of optimization theory and applications, 109(3):475, 2001.
Jialei Wang and Lin Xiao. Exploiting strong convexity from data with primal-dual first-order algorithms. In International Conference on Machine Learning, pages 3694-3702 PMLR, 2017.
Markus Weimer, Alexandros Karatzoglou, Quoc Le, and Alex Smola. Cofi rank-maximun margin matrix factorization for collaborative ranking. Advances in neural information processing systems, 20, 2007.
Markus Weimer, Alexandros Karatzoglou, and Alex Smola. Adaptive collaborative filtering. In Proceedings of the 2008 ACM conference Recommender systems, pages 275-282, 2008a.
Markus Weimer, Alexandros Karatzoglou, and Alex Smola. Improving maximum margin matrix factorization. Machine Learning, 72:263-276, 2008b.
Hongyi Wen, Xinyang Yi, Tiansheng Yao, Jiaxi Tang, Lichan Hong, and Ed H Chi. Distributionally. robust recommendations for improving worst-case user experience. In Proceedings of the ACM Web Conference 2022, pages 3606-3610, 2022.
Jason Weston, Samy Bengio, and Nicolas Usunier. Wsabie: scaling up large vocabulary image annotation. In Proceedings of the Twenty-Secon international joint conference on Artificial Intelligence- Volume Volume Three, pages 7764-2770, 2011.
13

---

Published as a conference paper at ICLR 2024
Huifu Xu and Dali Zhang. Smooth sample average approximation of stationary points in nonsmooth stochastic optimization and applications. Mathematical programming, 119:371-401, 2009.
Peng Xu, Jiyan Yang, Fred Roosta, Christopher Ré, and Michael W Mahoney. Sub-sampled newton methods with non-uniform sampling. Advances in Neural Information Processing Systems, 29, 2016.
Hsiang-Fu Fu Yu, Cho-Jui Hsieh, Si Si, and Inderjit S Dhillon. Parallel matrix factorization for recommender systems. Knowledge and Information Systems, 41:793-811 819, 2014.
Mofreh R Zaghloul and Ahmed Ali. Algorithm 916: computing the faddeyeva and voigt functions. ACM Transactions on Mathematical Software (TOMS), 38(2):1-22, 2012.
Yunhong Zhou, Dennis Wilkinson, Robert Schreiber, and Rong Pan. Large-scale parallel collaborative filtering for the netflix prize. In Algorithmic Aspects in Information and Management: 4th International Conference, AAIM 2008, Shanghai, China, June 23-25, 2008. Proceedings 4, pages 337-348. Springer, 2008
14

---

Published as a conference paper at ICLR 2024
Appendix
Table of Contents
Convex and separable upper bound pairwise loss
15
B Implicit Alternating Least Squares (iALS)
16
C Convolution-type smoothing
17
C.I Definition and derivatives
17
C.2 Interpretation
C.3 On the convolution type smoothing for CVaR
Properties of convolution-type smoothed check functions
Implementation details of SAFER:
20
D.1 Convolution-type smoothed quantile estimation
20
D.2 Re-weighted ALS
21
D.3 Instantiation of SAFER:
22
D.4 SAFER,+
23
E On Tikhonov regularization
24
F Alternating direction method multipliers (ADMM)
27
Details and additional results experiments
27
Models...
27
G.2 Evaluation protocol
28
G.3 Additional experiments
29
A CONVEX AND SEPARABLE UPPER BOUND OF PAIRWISE LOSS
We first derive the differentiable upper bound of the loss function in Eq. (1).
1 Σ L E HJl(i,)') fo(i,))),
M
JEV, PEV
max(0,fo(i.)) - fo(i.5) +1) JEVi''ev
M
Σ Σ
- +12
M
JEV
Note that used (a + 1)2 > 1 fxr
0} to obtain the third inequality.
15

---

Published as a conference paper at ICLR 2024
We next derive separable upper bound of the above loss follows:
E Σ [fa(i.3')
fo(i.5) 11
M JEV.
Σ Σ (1 (1-1a((,))) + 5o(e.519))
V
J'EV
[VI - (P(((1-)(( Σ 5o(1,))
M
JEV.
JEV
1
Σ JEV.
[1 1 fo(i.)))
H M
Σ fo(i.)) jEV
(16)
[
where μ 2 a hyperparameter. By substituting 8o = 4/lV| 2 1/111, we obtain our proposed
loss.? The obtained loss is convex w.r.t the predicted score fo(i. and separable w.r.t (t. and fo( if + j'. It also implies that, when using an MF model, fo(i.)) = (u,, vy), the loss is separable and block multi-convex w.r.f. u, Vy and Vyr for all (j.j') V22 such that j H
B IMPLICIT ALTERNATING LEAST SQUARES (IALS)
In this appendix, we provide brief description of iALS (Hu et al., 2008; Rendle al., 2022). The objective of iALS can be considered as a variant of ERM-MP, which is defined as follows:
min U.V
HALS(VUL V) + Q(U, V).
(17)
iE
where
GALS Vu, V.) = Σ v, + 3 Vu
(18)
JEV.
Observe that iALS's loss function is unnormalized for each user in contrast to the loss in Eq. (4). Using this, iALS solves the above optimization problem by alternatingly solving the block-wise subproblems as follows:
U(+1) = arg min
Σ iell
GALS (V Vi++ A 1
U
(19)
F
U
v(+1) arg min
Σ (VL (++1) V2) + || /1 | Veul
(20)
V
Each subproblem is row-wise separable, and we can obtain a closed-form update for each row U and V by solving the following linear systems,
Σ JEV.
+ Bu
+
+
u, Σ
(1)
(21)
JEV
JEV.
Σ iEL
(++1)
u
(++1)
+
30
(++)
(t- 1)
+ x
v, Σ u
(++))
(22)
EEU
REU,
where x 8 x = xx is the outer product operator. Various regularization strategies have been proposed for iALS (Rendle al., 2022; Zhou et al., 2008). Rendle et al. (2022) proposed the following weighting strategy,
(A.)., 1 X (DV/ 3,11)), =
(R4) -0141))
(23)
In this paper, we consider
as recently suggested by Rendle et al. (2022).
TThe second term on the RHS of the final inequality is known the implicit regularizer (Bayer et al., 2017)
and the gravity term (Krichene et al., 2018).
16

---

Published as a conference paper at ICLR 2024
C CONVOLUTION-TYPE SMOOTHING
C.I DEFINITION AND DERIVATIVES
We define the ramp function p (-) and the check function p-(-) for E (0, 1) as follows:
0, u,
0
and Pr(u) (T- 1)u,
u S 0 u >>>> 0
PI (1(u)=
(24)
TU:
Given a E (0. 1), we consider the convolution function (p1-a * kn) with the kernel density function kx(+) h- (-//b) with bandwidth h >> O. We can show that
(01-a (w) / pi Pi-al (n-a()) u)dv
kn u)du + (1 (1-a) l0a kn(u u)do
v. kn - uu)d + f. (re V- 1--kal u)dv
(t + u) kn (t) [ (t (t+ (t)
where an application of a change of variables yields the last equality. Then, we can obtain first derivative of the convolution function as follows:
V Vu(Pi-a * kn)(u) -a L ko(t)dt f kn (t)
= = + (1 (1- (1-1ca kn (a(()d))
= (1 a) -K(l-u).
(25)
Also, we can obtain the second derivative
+ = (n(--u).
(26)
C.2 INTERPRETATION
In this subsection, we offer a general understanding of convolution smoothing as a kernel density estimation for the underlying errors. Let y be continuous random variable having the density function fy- We consider the problem of predicting y by parametric model go(x) with parameters f and predictors Given some function p T [0, oo), consider population minimization problem,
min
- -(p(y-9((z)))
(27)
Applying the change of variables, can show that
Ejp(y ge(x))) = / ply (a--g(=))ff(()ee
/ plu)fulu;))du
(28)
where fu(v;0)= fa(u+ go())) can be considered as density of u(() y (x)).
Given finite samples {(r 3: with the sample size being n, we can non-parametrically estimate the density (u; by the kernel estimator, defined by
(v
kn (()),,
(29)
17

---

Published as a conference paper at ICLR 2024
with kernel density function k+(+) == k-/hh) with bandwidth h 0. Then, write the finite-sample counterpart of the objective function as
I p(u)fu(v;00)d
p(o)ka(v ((o))a(u-u(o0dda (()))dv
==1
Σ kn),((()))
(30)
where the last equality is due to the definition of convolution-type smoothing. The result above show that the convolution smoothing can be interpreted as kernel density estimation technique for approximating the distribution of underlying errors.
C.3 ON THE CONVOLUTION-TYE SMOOTHING FOR CVAR
Proposition C.I. If the kernel function kn is symmetric, the subproblem of is equivalent to the following,
min Ilera * (E(Vu,,,, Vi) Va) (2)-))] - }
where P1 (u) = a - 1(u < 0F)u.
Proof. We show that, if the kernel function kn is symmetric, the following equality holds,
+ - (p: * (y) - E) = + a
+
(0l- * ky) - E).
(31)
From the definition the convolution operator, we have
(pl ku) / max(0, An(u u)dv lm V. kn u)du,
(32)
which implies
= = + (or * (y
/e v. kn(w - E)
(33)
Here, we also have
(PI-a +k)axk))(11) = -a -a/ 00
v. kn (v u)dv + 1oe w. -n-k - w)de.
100
(34)
follows that
1 E + -p11-a kn) E)
-
a
y v. knkl fy--ya dv + 1 f v 1.-kal {y )dv
a Jo
(u v) ku(u Ev )) 1 l kn(u (v
dd
1 fo ku {))
(35)
where the second equality is due that 1 (n {y )d and the third equality is due the symmetric kernel. By combining Eq. (33) and Eq. (35), have
E + a
(pr (N 2 = + -----e + (pi-a kn) E
(36)
(pr kw à) y 1 ( kn)(g E)
+
(37)
a
which completes the proof.
18

---

Published as a conference paper at ICLR 2024
ERM-QE Decomposition.
The above observation immediately implies Proposition C.I,
+ E(er
(pr * (E((U, V.) -) [
((Vu,, V.) + al
E(p * (e)), E
a|u iEll
iEN
iell
(1-o(U.V.G)
ERM
Smoothed quantile estimation
(38)
This result is also interesting in the sense that the CVRR objective (I, into the objectives of ERM and smoothed quantile estimation. However, we not use the RHS for the update of U and V because the term of smoothed quantile estimation is not block convex because (u) is convex but decreasing in u < 0; See also Figure 2.
(U, V.g)
can be decomposed
C.4 PROPERTIES OF CONVOLUTION- TYPE SMOOTHED CHECK FUNCTIONS
We here discuss the properties of smoothed check functions assumed derive SAFER2
Property C.2. For any bandwidth 0 and any quantile level [0, 1). the smoothed check function (pi kn) satisfies the following properties:
(1) (p1 * k() is non-decreasing.
(2) (pi * kx) is closed if the loss function & and smoothed quantile & take finite values.
(3) (p1 * kx) is strictly convex if the kernel function kú, has the full support.
Proof. The first derivative of the smoothed check function is V. kx) =1 (n(-u), u), and therefore, it is non-negative based on the definition of the CDF K This implies (1). By satisfying the assumptions where (a) the function (p1-o kp. is continuous and (b) its domain dom is closed, (2) holds. Moreover, since the second derivative of the smoothed check function is the density function (i.e., (a (u) = (w)). , (u) > 0 holds for the domain of u, then (3) immediately follows.
* kn)
In (I), we confirm that the smoothed check function maintains the non-decreasing property. Conse- quently, the composite function (p1 kg)) preserves block convexity if is block convex. By using (1)-(2), the biconjugate of (p1 k2, is (p1 * kn) itself (Boyd and Vandenberghe, 2004), enabling the separable reformulation Eq. (11); note that the assumptions of finite & and & for (2) may not be stringent conditions in practice. (3) implies that the
method in the step may require full-support kernels, such as logistic and Gaussian kernels. However, we can still use the kernels with the support (-1, 1), such uniform and Epanechnikov kernels, by using sufficiently large bandwidth, which allows us to expand the support on (-h, ht; we further discuss the instantiation of SAFER, with such kernel functions in Appendix D.3.
19

---

Published as a conference paper at ICLR 2024
Algorithm 3: SAFER2 solver.
Vi E,, u, ((ū,(o//V)))), VjEV, Vy ~ (o0.(o/VV)) (o/Va))aa g-0
for farttil 1 to do -Ejvv
Gv + VTV
= VVV-
Compute Vi E U, 4 2/1.1 (, Ao Gyu,
for 1 to I. do
Uniformly draw Us Cu
(U..E))
+ -
7 arg g)(00)) ?
s.t. V I- (U, V, -
t { 7. d
end for
f 1 to [ do
- 1 KE - Li)
EyEn, Vov] + GB + aaupA1 Ia) " [ Eyen
u,
31 M
end
Gu t U diagMat(z)U for j 1 to [V/ do
At V, - (Eiez; M.
u + BrG + a/U -1 "
E'eu, 24 [V. U;
end
end
D IMPLEMENTATION DETAILS OF SAFER,
This appendix provides a detailed description of SAFER2, including its various instances for kernel functions. Furthermore, we shall detail a variant of SAFER, for a large embedding size, utilizing the subspace-based block coordinate descent recently introduced by Rendle al. (2021).
Alternating optimization.
SAFER, updates each block cyclically as follows:
((++)) arg min
Σ (pr * kin (e) VV) - E)
{ j.
(39)
a/d1
VEU
2 (1+1) arg max
(((V(0)) V:) - -((+1)) - (pr kin) (#:) }
(40)
iEll
U(++1) arg min U
1
E( ((V()), u; D2) (I+1)
+ + 1/ U 1
(41)
iEll
((+1) ((Vu((t)) v.) + 1 2
v(+1) = arg min
(42)
a
V
Below, present efficient implementation of each step. The overall algorithm is described in Algorithm 3, including the update formulae for U and V.
D.1 CONVOLUTION-TVER TYPE SMOOTHED QUANTILE ESTIMATION
Wowton-Raphson algorithm. The subproblem of & for the general kernel density function b cannot be solved analytically. Hence, resort numerical solution using the efficient Newton-
20

---

Published as a conference paper at ICLR 2024
Raphson (NR) method. At the /-th iteration in the (t 1)-th update, estimate as follows:
Pe v()
< = (((+1) 7 aet11, where d(+1)
(43)
v 1- V
Here, the first and second derivatives of a) (1-a)-Cts-CVa) can be evaluated as follows:
VeW1-o a/
Σ [11 a) Ka(5 ((V,, Vi)l.
(44)
iEll
(U all
(n ((Vu, V.)).
(45)
iE
Pre-computing the loss ((Vu; V:) for each user can be helpful in reducing the computational burden when the gradient and Hessian can be computed exactly. Furthermore, an effective initialization of is crucial, and we set s to the previous estimate, i.e. < =
Backtracking line search. The value of 7 plays vital role in obtaining an accurate solution for Ç. However, constant y often fails provide efficient results. One way to overeome this is by employing the widely-used backtracking line search to determine & adaptively. It aims to find the maximum y E (1/2, 1/4,...) such that the Armijo (sufficient decrease) condition (Armijo, 1966) satisfied. We can express this as:
"" max 7E/1/2.1/4....
7
(46)
s.t. I
cyd Velevaaaaaa V, 7
where e >> 0 is the error tolerance parameter, often set to small value such as 10-4, Performing each iteration of backtracking line search for evaluating 11-a V. -5-- 7 d) takes the cost of ((1()). Therefore, this step demands cost of (.) for each epoch.
Sub-sampled algorithm. Although the computation of the direction d can be done in parallel for users, it may be costly for many users, particularly when using backtracking line search with the cost of L). We can introduce the sub-sampled Newton-Raphson method (Roosta- Khorasani and Mahoney, 2019) to reduce this cost by approximating the gradient and Hessian based on the uniformly sub-sampled users. Let [24, be the sub-sample size of users, which is smaller than the original users size (2]). Chernozhukov and Fernández- Val (2005) obtain the large-sample properties of the quantile regression estimator based on sub-sampling, under the setting where [44 [M/l/M - 0 [4] - Their results suggest that the sub-sample estimator & of the true value
and [24,] Eo satisfies that & to + (p/ (261). (1/v/)))) That is, the sub-sample estimator converges to the true parameter at rate of 1/ / [26)-- In practice, the user size often larger than ten million or more, and the sub-sample size [ 100,000 ensures that the estimation error asymptotically vanishes as 1/ / 0.00316.
D.2 RE-WEIGHTED ALS
The update of and can reformulated by using primal-dual splitting as in Eq. (11). Here, we focus on the update of primal variables, i.e., U and V, since we have already described the step in Section 3.3,
U step. Given z, the optimization problem of U and V forms re-weighted ERM. Owing to separability, the update of U.
U(++1) = arg min
1
Σ 12
(+1))
((V()) Vi) + IIA!
(47)
alu
U
iel
21

---

Published as a conference paper at ICLR 2024
can be solved in parallel with respect to each i, as follows:
((+1)
V &l (VV Vi) + ( u =
aM
((+1)
2(1+11 2f
(++1)
30 Σ
taRMIAST)
JE
JEV
JEV;
(b)
(48)
The updated u, can be obtained by solving the linear system above. Since the user-independent
Gram matrix indicated by (b) has been pre-computed (at cost of O(V]dd)) the computational cost
of updating each user's u, involves computing the user-dependent partial Hessian indicated by (a)
in ((VV.dd)) and then solving a linear system of d x d in Odd*). Since Eieu [V] =S|, the total computational cost is thus O(S d [4)
V step. Analogously, the update of V.
v(++) arg min
2
(++1)
((Vu
(++1) V.)l 1 1
(49)
[|
can be solved as follows:
Vv; a|
Σ
..V,, Vi) + 1 =0
VEU
(f+1)
(1+1)
((+1)
(I+1)
30
Σ
(1+1)
((+1)
e u +all11l01
((+1)
2j Σ SEU,
((+1)
Vi
[
Seu
(50)
We can reduce the computational cost of this step as in the U step by caching the weighted Gram matrix indicated by (c), which costs O((())2)) when the loss for each user is pre-computed. The computational cost for updating V is O( + Vld").
D.3 INSTANTIATION OF SAFER,
instance SAFER, is determined by the choice of the kernel density function kn- We shall describe the implementation of SAFER, with some popular kernels as it would be helpful for reproducing our method; the implementation described here available https: /ithhb.. safer22 orrrrorreeeeeeeeeeeeeeeeeeeeet recommender.
Gaussian kernel. For the Gaussian kernel, the kernel density kn and its CDF K₁, can computed as follows:
u²
kn)
exp
(51)
/2Th
2h.
1
K2 (u)
+ erf
erfc
(52)
2
/2h/
/2h
where crf(u) 2 / S exp(-v2)ddu and erfo(u) crf(u) are the error and complementary error functions, respectively.
These complex functions generally implemented as special functions and be computed very efficiently (Abraroy and Quine, 2011: Gautschi, 1970; Poppe and Wijers, 1990; Zaghloul and Ali, 2012).
The smoothed check function (PI-o * kn.) is then obtained as h [h ku)
(pi-a ku) which is used for backtracking line-search in the step.
2- kn(-w)) + ((1 a) 0.5)u,
(53)
22

---

Published as a conference paper at ICLR 2024
Epanechnikov kernel.
We also describe the implementation for the Epanechnikov kernel (Epanech- nikov, 1969). Epanechnikov kernel density kn and its CDF K₁, can be computed as follows:
"
27 1 1(//// s 1
(u) = 4h
(54)
h
0,
< -1 u/h [-1, u/h >
1
[u[3h2 - =²² +2h101,
(55)
1h
1,
The smoothed check function (Pi-a * kn)) is then obtained as
h
u
1/ (
(01- * ku)
2
4
T.
+ 8
11 H/wu/b/ /l 1
/ 1 1} + ((1 - a) 0.5)u.
(56)
Note that the support of (xx is on (-h, ht), and thus we can ensure the strict convexity (i.e., positive Hessian) of (pi-a kn)) in Newton-Raphson step.
D.4 SAFER
The SAFER, algorithm experiences the quadratic/cubic runtime dependency on the embedding size d. This problem has been tackled by various studies (Bayer et al., 2017; He et al., 2016; Pilászy al., 2010), and some studies recently reported that the large dimension is very important to improve ranking quality of iALS (Ohsaka and Togashi, 2023; Rendle et al., 2021) Considering this, we propose an extension of SAFER, for a large embedding size by using the recent subspace-basd basedd block coordinate descent of Rendle et al. (2021), which is a simple yet effective approach, enabling efficient utilization of optimized vector processing units.
Subspace-based block coordinate descent. To overcome the above problem, iALS+ considers the subvector of user/item embedding as a block (Rendle et al., 2021) and optimizes the subvector by a Neewoon-Raphsee method. We here apply this approach to SAFER, Let S {1, {1...., d) be a vector of indices and u,, be the subvector of u, corresponding to æ. Then, the first and second
A53 |u: in Eq. (11) with respect
derivatives of the objective [u(,; V) = (z:/a[(() e(Vu;, V3) + to ur are:
Vu, LU V)
Σ CvJJus u, 1) + 2,30
Σ V
u, + a u
(57)
Vi
jEV;
JEV
V2 Li V)
V J,T 2100 Vj
Σ
i + 1
(58)
JeV.
JEV
Note that omitted the constant factor (a(21))1 for brevity. We cre-compute the partial Gram matrices Ejev V in ((V][T]d)) Ejev V,, in ((V][T ), and the prediction vjJ u, for (i.j) E S in ((S]d). Then, the computational cost of the gradient and Hessian for all users is O(S|| ||| + [2]]] Rld + (²) We subsequently update the subvector u by a Newton-Raphson step,
(++1)
-V Lu v(U) -V n.. L.(
(0) V(O).
(59)
This can be computed by solving linear system of size .| | in O(1 | ()) time. The computational
cost for updating all user subvectors is O(S ग. + |V d + S|| 2 + 13))
Similarly, denoting (U, v) = (a - EiEula ((Vu, WJ] + 2
| 12, we can obtain the
update of an item embedding as follows:
(1+1)
j. - (V? (V2 L, V.. L.U,V)
(60)
23

---

Published as a conference paper at ICLR 2024
Algorithm 4: SAFER,+ solver.
Vi EE N(O, VjEV, Vy ~ (o/VV)L)
g-0
for torttll 1 to do -jvvv
Gv VVV- =
VTV
Compute Vi E U, 4 2/1.1 (1 Ag Gyu
for k to I. do
Uniformly draw Us Cu
(U..E)
+ -
γ arg g)(00)) ?
s.t. V I- (U, V, - 7d)
t 7. d
end
for f 1 to [ do + (E (i)
end
for π P do
G for
V G EjEv V i
1 to [1] do
, u, (Ve u L,(,, V))- uin Lu V)
end
Ggl + Σieuh U u: G + ΣiEu 24 U1 U.
for j f to | do
V f V (V V3. L,(U, - VV L,
end
end
end
where
Vv,, (), σ
24
(u V 1)v + +.30
Σ 2i u v, + a V (61)
[4
iE44
ViEl
V2 - L,(U,) x
Ug, 1, + +.0
I
(62)
[V
The computational cost for updating all item subvectors is O(SS| + [ |
V|[ (3)
Computational complexity. We follow the iteration scheme suggested by Rendle et al. (2021), which cyclically updates the subspace of user and item sides for each subset of indices. As result, we obtain the following computational cost.
o
SS|d 4 + (nu) + []] + | + |m|) + |
+ +|
(63)
O(S (|[S|T// + (I + (d + )).
(64)
The algorithm shown in Algorithm 4.
E ON TIKHONOV REGULARIZATTION
In this appendix, we develop a regularization strategy for SAFER,, which allows us to control the numerical stability of subproblems for users and items. Since setting appropriate regularization
24

---

Published as a conference paper at ICLR 2024
weight for every user/item is impractical in large-scale settings, derive single hyperparameter that controls the regularization weights for all users and items.
For matrix A E _ddd the condition number x: x(A) is defined follows:
K( (A) = A || ||- ||
(65)
The condition number of a matrix A characterises the numerical stability of linear system Ax = b; this problem with respect is numerically unstable when x(A) is large. We consider normal A, and then the condition number can be computed as follows:
Amax (A)l
(A)
(66)
Ammin (A)]
where Amax (A) and Amin(A) (A) are maximal and minimal eigenvalues A, respectively. In our case, want keep small the condition number of each Hessian matrix with respect to the rows of U and V. For the i-th row of U, the Hessian matrix H, is as follows:
H:
= Zi
Σ JEV
V
+ 2, V V + I
(67)
which implies
Amia H.)) a/Z
Amax (H.) Amax
HH a ( a/
(68)
To ensure a small value of <(H.), ), we can adjust the regularization weight A > 0. However, the maximal eigenvalue Amax (H,) changes at each update step in the alternating optimization, but computing the dominant eigenvalue is a costly process (Mises and Pollaczek- Geiriiger, 1929). Therefore, propose simple regularization strategy control the condition number of each linear system solved with a constant hyperparameter. Assuming v > O is the upper bound of the squared norm of the user and item embeddings throughout the optimization, i.e., v 2 1|v for all j E V and
2 | 22 for all E U, we have the following upper bound for Amax (H; a [2 A I):
Amax (H: a(( <
Vjv
2,80 V
Vi
jEV:
1vy 11 + 2: VI
[ JEV:
$ 2, + 2,3p 1/P
S v(1 Vo)
In the last inequality, we used z; S 1. Therefore, by setting
X
(
(1 3ol
(69)
a
with hyperparameter À > O, we can ensure the following incquality,
Amax H. al I) + Amin (H,)
K (H))
Amin (H;)
< < +1.
(70)
The advantage of this reparametrization is that we may be able to bound the condition number of each user's linear system from above by (w/X) + which independent of ax, [24], [V. and This allows to flexibly control the regularization intensity through tuning À while ensuring that the subproblems of each user are conditioned to the same extent. Note that too- value of À leads to poor model training while the condition number will be close one, and therefore, still need to tune X.
Analogously, we can derive the regularization weight for the j-th row of V.
H =
2;
n + UT diagMat(2)U + {
(71)
iEM
25

---

Published as a conference paper at ICLR 2024
and
Amax (H - a
Σ
2
n
+ U diagMat(z)U U
[il
1E
F
2
|u
12 + 3o
2/
u 112
[
+ 30
iEll,
iEU
+ 30
Σ
iEM
VEU
In contrast to the case of user embeddings, we applied S 1 for the first term of the last inequality. This is because bounding Lieu 2; S [4] rather loose, and the weighting strategy based on this bound will lead to the over-regulariiation of item embeddings. To avoid this, we introduce the following property of convolution-ty type smoothed quantile.
Proposition E.1. Suppose that samples {li,...,fy )n) of losses and its convolution-type smoothed quantile En. == arg ming \ (p1-a kn) E). Then, n dual variables {z1,..., Zn} satisfy
Σ 2
(72)
where 2; K lz).
Proof. The smoothed quantile En satisfies the first-order optimality condition,
Pi-a ke) En) =0
[[1 a) -Kn(En - -4)]
Knlen Ec] = =a-
which immediately completes the proof.
From this result, we can substitute Eier 2, in the upper bound of Xmax (H, a((()1 with a/4 and then obtain
Amax (H aal((1)) Su
Σ iEU,
Ddl + 8ocil21
(73)
By setting
+ B
(74)
iELl
we can ensure the following inequality
k H = Amax (H, a/ (2
Amin (H))
Amin (Hy)
v( H ee, MTT + Soa(2()
+1
M iie; D. + (1)
+ MEER, [&t + 3ua(U1)
A +1
26

---

Published as a conference paper at ICLR 2024
On the regularization strategy of iALS. Tikhonov regularization has been widely adopted for MF models with the ALS solver (Rendle et al., 2022: Zhou et al., 2008). In particular, the recent technique in Eq. (23) (proposed by Rendle et al. (2022)) can be obtained by similar derivation. Namely, consider the Hessian matrix and regularization weight for the i-th user,
H. Σ v.v +VVVV+++ + JEN
(75)
A ([V)) + Bol
(76)
then we have
Amax
H - S Σ
Vov + V
V Vil
JEN
2 Σ
1 + V
JEV
$ < + B V1).
which implies s:HH) v/ + 1. The result for each item is analogous, and we therefore omit the derivation.
F ALTERNATING DIRECTION METHOD OF MULTIPLIERS (ADMM)
As discussed in Section 3, the smoothed check function (01 kv)) non-linear and leads to the coupling between the rows of V as follows:
min U.V.
E+ 0.41
E(er * kin) (E(Vu,, V.) - -) 1 VEU
1 11A
To decouple the rows of one can consider the use of the alternating direction method of multipliers (ADMM) (Boyd et al., 2011; Steck et al., 2020: Togashi and Abe, 2022) by introducing auxiliary variables E RR/4, which leads to the following constrained optimization.
1
(pr kn) (w -) + 1 A +
min U.V.
E
(77)
a/
s.t. yi (Vu,, Vi), Vi E u.
The augmented Lagrangian in scaled form is defined as follows:
(U, V.E.y, w) = ax/21 Σ (pr kn) (ye E) 1 A U + 1
1
1
P (wy ((Vu, V.) y.)
(78)
where E MM is the dual variables (i.e., the Lagrange multipliers). Observe that, because of the quadratic penalty term of ADMM, the rows V are still coupling in the objective. One can avoid this by using another reformulation by introducing Iim-einenso auxiliary variable Vu; for each user, which requires prohibitively large dual variables of size [|[[)
G DETAILS AND ADDITIONAL RESULTS OF EXPERIMENTS
This appendix provides detailed descriptions of the experimental settings and additional results, which are omitted for the strict space limitation.
G.1 MODELS
iALS. Our implementation of iALS is based the reference software publicly provided by Rendle et al. (2022). This iALS implementation is reported to be competitive with state-of-the-ar methods on ML-20M and MSD datasets. We used their proposed regularization strategy with = 1.0 as
suggested by Rendle et al. (2021). We also followed their implementation of iALS+ and set v 1.0 for all settings as in iALSS.
27

---

Published as a conference paper at ICLR 2024
Table 3: Ranges of hyperparameters.
Models
Hyperparameters
iALS
5g E 1,), E 5e-4,0.1
Mult-VAE
T [1e-4, le 2]. [4 [50, 100, 300, 400, 500), e [0.1, 1.0]
2
ERM-MF
Bu E 1c- 2 [1c-4,
CVaR-MF
Bg E (lc-3,0.1), E [le-4,1e-2] 0 ( 0.51
SAFER,
Bo E [1e-4, A (le-4,1e-2,a=033, 0.3, [0.1, (0.1,1.0) 1.0), RMI/// 0.1, -= =
ERM-MF and SAFER2. We implemented ERM-MF and SAFER, (SAFER2*t) in the same code- base as iALS. For ERM-MF and SAFER,, we used our proposed Tikhonov regularization as we found it is generally effective in terms of the final quality and hyperparameter sensitivity.
CVaR-MF. The implementation of CVaR-MF is also provided in our software. As CVaR-MF often takes much longer time converge, we tune constant learning rate found that applying the gradient descent to the & step is quite unstable and makes it hard to obtain acceptable performance. Therefore, we exactly compute the (1 a) -quantie of the users' loss as in each step. To finely observe the difference of the solvers, we used our proposed regularization also for CVaR-MF, which empirically leads to good performance.
for U and V. We also
Mult-VAR We implemented the method based on the public codebase." We tune learning rate 7 > 0, batch size (4), and annealing parameter 0.
Prediction for users. As follow the strong generalization setting, each MF-based model must produce predictions for new users who are not in the training split and thus do not have the trained embeddings (e.g., uz). To this end, we follow previous studies (e.g., Rendle et al. (2022)), where each model solves an independent convex problem for each user by leveraging the 80% of the user's interactions. In iALS, we can obtain the embedding u, of new user with VA by solving the following problem:
(AALS(VU V.)
ug = arg min uERd
1
(79)
where X X(D2) + 3o/11).
In ERM-MF, CVaR-MF, and SAFER2, the problem can be expressed as follows:
X
u; arg min uERd
al (Vu, V.)
| 2
(80)
where set x (A/a + (o]V),, and V is the trained item embedding matrix; fix it in the prediction phase. This is standard least-square problem and hence easy to solve by computing the analytical solution. Note that applying this prediction procedure CVaR-MF is reasonable because the subproblems for users are completely independent in this step.
the gradient-based
Hyperparameter tuning. In the experiments in Section we tuned all models by using the validation split for each dataset. The range of each hyperparameter is presented in Table 3. We tuned all the hyperparameters by performing grid search for ML-IM. For ML-20M and MSD, we tuned all the parameters manually to reduce the experimental burden. The number of epochs T is set to 20 for iALS, ERM-MF, and SAFER; in the hyperparameter search and set to for training the final models. For CVaR-MF, we set T 500 for validation and T 1,000 for testing.
G.2 EVALUATION PROTOCOL
Datasets and pre-processing protocol. We employed a standard pre-processing protocol for the datasets (Liang al., 2018; Rendle al., 2022; Steck, 2019b; Weimer et al., 2007). The
httpps://i thwb. httpss/////ttto.ooooooooooeeeeee pyytor
28

---

Published as a conference paper at ICLR 2024
implementation of the pre-processing protocol is based on Liang et al. (2018). As described in Section 5, we divided the users into three subsets: the training subset (i.e. {V; heu) contains 80% of the users, and the remaining users split into two holdout subsets for validation and testing purposes: For each validation and testing subset of MI.IM, MIL-00M, and MSD, the number of users evaluated is 1,000, 10,000, and 50,000, respectively.
Evaluation measures. In our experiments, we use recall at K (R@ K) and normalized discounted cumulative gain K (nDCG@ K) as measures of ranking quality, Let V/ be the held-out items pertaining to user i, and n; (k) E } be the k-th item on the ranked list evaluated for user i. The computation of R and DCG@ follow:
K [lm EV).
RaK(m, V2) min(K, [V')) k=1
(81)
K
1 (k), V
DCG@K(T, = Σ k=1
(82)
loga(k +)
DDGG@ K is defined as nDCG@K(T;; V) = ( where π is an ideal ranking for user i.
G.3 ADDITIONAL EXPERIMENTS
Effect of the sub-sampled Newton RewtonnRaphson method. Figure shows the effect of the number of sub-samples [66] for each sub-sampled NR iteration the step. Each curve was obtained by varying [24]] with the best hyperparameter setting for [44l///] 0.1. We can observe that (1) the primal variables (i.e., and V) converge even with a small [461; (2) the dual variables (i.e., z) fluctuate with small (4, values; and (3) the final ranking qualities almost identical for different [44] values. It suggests that the sampled NR method effective in practice as alleviates the computational cost of the step, which is the only additional cost from iALS.
We also report the effect of the number of iterations trend in Figure 7: Small values of L and U4, lead to the fluctuation of whereas the convergence is maintained in most cases. The final quality slightly deteriorates when both L and (24 are small, particularly in terms of the semi-worst-case performance (i.e., a 0.3).
the & step in Figure 8. There is a similar
[6//MM|
[6/// 24| 0.05
R4////4
[ 1.0
|
VV-1 ViIle
g
Resicual norm
Residual norm
10"
101
10F
10
Residual
10-
109
10:
10
40
10
20
30 Training epochs
Training epochs
Training epachs
1.0
0.3
0.855
0.940
R01100
0.650
p88100
0.386
R00100
0.15
0.845
0.14
0.325
0.640
10
20
30
50
0.320
10
20
30
40
0.15
10
30
40
50
Training opochs
Training epachs
Training epochs
Figure 7: Convergence profile of SAFER, with different sampling ratio [61//4 on ML-20M.
Choice of kernels. Various symmetric kernels can be used to instantiate SAFER, discussed in Appendix D.3. To observe the effect of choosing kernel functions, we compare SAFER, with the Gaussian kernel, examined in Section 5, and with the Epanechnikov kernel. Figure 9 demonstrates the distribution relative ranking performance of each method compared to iALS through nested cross-validation in Section 5. The SAFER, instance with the Gaussian kernel performs slightly
29

---

Published as a conference paper at ICLR 2024
[6 0.1,
1.0,
IV [66///
0.1,
Rab//14 1.0,
1000 U'lle
110
Residbal norm
101
Rosiual norm
10"
Residual
102
100
10'
10
20
40
80
0
20
40
80
0
20
40
80
Training epochs
Training epuchs
Training epuchs
0.635
1.0
0.3
0.1
0.340
R8150
0.650
P8130
0.335
R01000
0.15
0.645
0.14
0.325
0.840
40
0.320
0:13
Training epachs
Training epacts
Training epacks
Figure 8: Convergence profile of SAFER, with different number of NR iterations L on ML-20M.
Relative performatce vs 1ALS
RC20 1.0)
RE20
0.3)
Reaaiive prrformttce 95 1ALS
R050
1.0)
Relative performance 15 IALS
R850
0.3)
1.00
Relative performaace 95
5.10
1.00
7.05
1.0
1.z2
1.00
1.00
1.00
3.66
BRM-MP
SAFER, Fjm
SAFER Ganes
BRM-MP
SAFER Rjn.
SAFER
ERA-MP
SAPER Rje.
SAFER Classe
ER3:MP
SAFER, Rja
SAFER
Clanss
Glrss
Figure 9: Distributions of relative performance vs iALS on ML-1M.
better than the one with the Epanechnikov kernel, but the difference between them is not substantial. However, the Gaussian kemeel results in stability in the & step due to its full support property, making it easier to tune the bandwidth parameter h.
Inexact linear solvers for large embedding size. Because SAFER, has quadratic/cubie runtime dependency on the embedding size d, we proposed a variant of SAFER,, i.e., SAFER,*+, in Ap- pendix D.4. Here, we examine the computational efficiency of SAFER To establish a baseline method for comparison, we consider another variant of SAFER,, called SAFER;-CG, that uses the conjugate gradient (CG) method to solve d linear systems. For SAFER;-CG, we used the CG implementation in the Figen library"; the maximum number of iterations was set to five, and the error tolerance was set to 1e-4. We used the same hyperparameters for models, which achieved the best ranking quality with the exact linear solver and d 256 for the validation split of ML-20M. In Figure 10, we present the convergence speed SAFER,* on ML-20M for different values of d and |t|. For comparison, SAFER,'s results, except for the case of d = 4,096, where SAFER; did not finish in a practical time. Our results show that SAFER,+ achieves comparable ranking quality compared to SAFER, for both average and semi-worst-case scenarios. Furthermore, although the convergence speed in terms of training epochs similar for both SAFER; and SAFER,*+, SAFER exhibits substantially superior computational performance. When d is small (c.g.. d = 256 in the second row of the figure), SAFER;-CC (light green line) outperforms SAFER,+ terms of wall time. However, the performance gain of SAFER,* increases for larger values of d, such as d = 1,024, 2,048, or 4,096, highlighting the scalability of SAFER,* with respect to the embedding size.
also display the red dashed curve that corresponds to the original
Phttps://eighe www.aaatntttinnmmemm oy.org
30

---

Published as a conference paper at ICLR 2024
[F]
10
-1/=32
64
[r] - 125
SAFED
SAFER,-CO OG
0.30
0.50
Tioking easts
Traiing spocks
Wallirs (ses)
Wall siw (uc)
236,
350
0:21
0.51
0.30
Troking epache
Tox spocks 312 00
Welskee (ee)
Well siee (te)
di 518 16
510 11
a
0:00
0.60
Triking
Trasig
d 1024 u.
d 1104 a 1tS
1024
- 1104
0.00
0.24
0.00
(IN)
(sc)
df 0118 o 1n
d 2018 0.
di 2045 m- 1n
of c0us, U.S.
0.64
80000
0900
8.56
0.90
024
10
Triiing spods
spock
Wall tire (eec)
Wall sime (wc)
J 32%
d UCG 1LS
Jl 4235
- 4186, 0.5
0.44
0.64
0.02
0.00
0.50
0.50
127
17
107
Traiiig saoda
skg spocts
Wylikce (IN)
Will Sime (wc)
Figure 10: Effect of embedding and subspace block sizes on convergence speed of SAFER
31


Reviews:
[{"review_id": "TBJ7vb0tO6", "reviewer": "Reviewer_nPPr", "summary": "The paper proposes SAFER2 for excellent tail performance. \\\nThe paper modified CVaR and devises CVaR-MF and SAFER2.", "questions": "Please refer to Weaknesses.", "limitations": "", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "1. Well-written and organized.\n2. Solid derivation and equations.", "weaknesses": "[Major cencerns]\n1. Motivation\n- As far as I understand, the main motivation of this paper is tail performance. (in Abstract)\n- However, I found that the main contribution of this paper focuses on parallel computation and Scalability. (in the manuscript.)\n- What if we just put some high weight on the tail users? or What if we just compute gradient only with bottom-$\\alpha$% losses?\n\n2. Experiment\n- The main experimental result focuses on that SAFER2 shows comparable results when $\\alpha=1$ and superior results when $\\alpha=0.3$.\n- When $\\alpha=1$, there are many recent methods outperforming Multi-VAE. e.g., LightGCN, SGCL. and in my opinion, the pair-wise and list-wise loss functions outperform the point-wise loss.\n- When $\\alpha=0.3$, there are many recent methods enhancing the performance of tail users.\n- Did you train two SAFER2 models with $\\alpha=1$ and $\\alpha=0.3$? If so, that is somewhat unfair for the baselines (e.g., VAE-CF).\n\n[Minor cencerns]\n1. Equation numbers are needed for each equation when finding the location of the equation by ordering.", "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": "", "extra_scores": {}}, {"review_id": "23TKL5wIxQ", "reviewer": "Reviewer_RZsQ", "summary": "Authors propose a method for recommender systems via Matrix Factorization, that is having a special type of loss-function that is connected to conditional value at risk (CVaR in Finance). \nHowever, CVaR is \"hard\" for optimization, author re-use Rockafellar and Uryasev (2000) reformulation of CVAR to block multi-convex  loss function, over which CVaR-MF optimzation can be done. Still this formulation is not scalable, since max(0,x) part prevents separability for items. Author propose to use Molifier kernels to get strict convexity. \nFinally, they wrap the steps in Smoothing Approach For Efficient Risk-averse Recommender (SAFER2).", "questions": "Can you make comparison to the existing model:\ne.g. Wen, Hongyi, et al. \"Distributionally-robust Recommendations for Improving Worst-case User Experience.\" Proceedings of the ACM Web Conference 2022. 2022.\n\nIn your paper, major contributions are done by forcing re-formulation of CVaR to be scalable. Can you elaborate, what is the problem of controlling the \"tail-risk\" with distribution-wise losses like divergences (KL), or Wasserstein metrics?", "limitations": "", "rating": 8, "confidence": 3, "soundness": 4, "presentation": 3, "contribution": 3, "strengths": "Going beyond standard ERM objectives, nice mathematical re-formulation of tail risk CVaR and combining it with Matrix Factorization for Recommender Systems.", "weaknesses": "No clear weakness, except more baselines from Robust Recommender Systems field.", "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": "", "extra_scores": {}}, {"review_id": "5Avu5gOaRy", "reviewer": "Reviewer_QW8m", "summary": "The authors target tail performance in recommendation systems- a crucial but often overlooked scenario in the recommendation domain. They propose a collaborative filtering-based training approach for minimizing the conditional value at risk (CVaR), which represents the average risk for tail users. This is in contrast to conventional ERM-based methods, which focus on the overall average performance. The authors first show that a reformulation of the ERM objective using matrix factorization includes a non-smooth and non-linear term, which hampers the objective’s separability, making it difficult to parallelize and scale the approach. The authors circumvent this issue by applying convolution-type smoothing in the objective. They then use convolution-type smoothed quantile estimation in their proposed method, “Smoothing Approach for Efficient Risk-averse Recommender (SAFER2). They evaluate their model in the semi-worst case and average-case scenarios and find that, in general, their approach improves tail performance while maintaining average case performance.", "questions": "No particular questions.", "limitations": "", "rating": 8, "confidence": 2, "soundness": 4, "presentation": 3, "contribution": 4, "strengths": "1. The authors propose a novel and efficient approach for making CVaR minimization tractable, allowing them to target tail performance in their training objective.\n2. Their approach improves tail performance without significantly impacting average performance.\n3. Rigorous mathematical proofs are provided for justifying each step of their algorithm.\n4. The authors compare the performance of their model with several baselines and state of the art methods.\n5. The experiments also demonstrate that the proposed method is robust and computationally efficient.", "weaknesses": "Some concepts may need to be explained in more detail in the main paper (e.g.- the reformulation of the CVaR minimization objective, smoothing techniques for quantile regression, NR algorithm etc.)", "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": "", "extra_scores": {}}]

Author rebuttal:

Decision: Accept (poster)