Scalable Learning and MAP Inference for Nonsymmetric Determinantal Point Processes
Authors: Mike Gartrell
Insu Han
Elvis Dohmatob
Jennifer Gillenwater
Victor-Emmanuel Brunel
Venue: ICLR 2021 Oral (2021)
Keywords: determinantal point processes
unsupervised learning
representation learning
submodular optimization

Abstract: Determinantal point processes (DPPs) have attracted significant attention in machine learning for their ability to model subsets drawn from a large item collection. Recent work shows that nonsymmetric DPP (NDPP) kernels have significant advantages over symmetric kernels in terms of modeling power and predictive performance. However, for an item collection of size $M$, existing NDPP learning and inference algorithms require memory quadratic in $M$ and runtime cubic (for learning) or quadratic (for inference) in $M$, making them impractical for many typical subset selection tasks. In this work, we develop a learning algorithm with space and time requirements linear in $M$ by introducing a new NDPP kernel decomposition. We also derive a linear-complexity NDPP maximum a posteriori (MAP) inference algorithm that applies not only to our new kernel but also to that of prior work. Through evaluation on real-world datasets, we show that our algorithms scale significantly better, and can match the predictive performance of prior work.

Published as a conference paper at ICLR 2021
SCALABLE LEARNING AND MAP INFERENCE FOR NONSYMMETRIC DETERMINANTAL POINT PROCESSES
Mike Gartrell
Insu Han KAIST
Elvis Dohmatob Critco Al Lab
Critco AI Lab
m.gartrell @criteo.com
insu. haaukaa@kai ac.
e. e.dohmatoon tobecri teo. com
Jennifer Gillenwater
Victor-Emmanuel Brunel
Google Research jengi Sangiggeogge..o. com
ENSAE ParisTech
victor. comanuel..bureeeeeeee
ABSTRACT
Determinantal point processes (DPPs) have attracted significant attention in ma- chine learning for their ability to model subsets drawn from a large item collection. Recent work shows that nonsymmetric DPP (NDPP) kernels have significant advantages over symmetric kernels in terms of modeling power and predictive performance. However, for an item collection of size M, existing NDPP learning and inference algorithms require memory quadratic in M and runtime cubic (for learning) quadratic (for inference) in M. making them impractical for many typical subset selection tasks. In this work, we develop a learning algorithm with space and time requirements linear in M by introducing composition. We also derive a linear-complexity NDPP maximum posteriori (MAP) inference algorithm that applies not only to our new kernel but also that of prior work. Through evaluation on real-world datasets, we show that our algorithms scale significantly better, and can match the predictive performance of prior work.
NDPP kernel de-
INTRODUCTION
Determinantal point processes (DPPs) have proven useful for numerous machine learning tasks. For example, recent uses include summarization (Sharghi et al., 2018), recommender systems (Wilhelm et al., 2018), neural network compression (Mariet & Sra, 2016), kernel approximation (Li et al., 2016), multi-modal output generation (Elfeki et al., 2019), and batch selection, both for stochastic optimization (Zhang et al., 2017) and for active learning (B1yık et al., 2019). For subset selection problems where the ground set of items to select from has cardinality M, the typical DPP is parameterized by an M x M kernel matrix. Most prior work has been concerned with symmetric DPPs, where the kernel must equal its transpose. However, recent work has considered the more general class of nonsymmetric DPPs (NDPPs) and shown that these have additional useful modeling power (Brunel, 2018; Gartrell al., 2019). In particular, unlike symmetric DPPs, which can only model negative correlations between items, NDPPs allow modeling of positive correlations, where the presence of item in the selected set increases the probability that some other item will also be selected. There are many intuitive examples of how positive correlations can be of practical importance. For example, consider a product recommendation task for a retail website, where camera is found in user's shopping cart, and the goal is to display several other items that might be purchased. Relative to an empty cart, the presence the camera probably increases the probability of buying accessory like tripod.
Although NDPPs can theoretically model such behavior, the existing approach for NDPP learning and inference (Gartrell et al., 2019) is often impractical in terms of both storage and runtime requirements. These algorithms require memory quadratic in M and time quadratic (for inference) or cubic (for learning) in M; for the not-unusual M of 1 million, this requires storing STB-size objects in memory, with runtime millions or billions of times slower than that of a linear-complexity method.
In this work, make the following contributions:

---

Published as a conference paper at ICLR 2021
Learning: We propose new decomposition of the NDPP kernel which reduces the storage and run- time requirements of learning and inference to linear in M. Fortuitously, the modified decomposition retains all of the previous decomposition's 's modeling power, as it covers the same part of the NDPP kernel space. The algebraic manipulations we apply to get linear complexity for this decomposition cannot be applied to prior work, meaning that our new decomposition is crucial for scalability.
Inference: After learning, prior NDPP work applies a DPP conditioning algorithm to do subset expansion (Gartrell al., 2019), with quadratic runtime in M. However, prior work does not examine the general problem of MAP inference for NDPPs, i.e., solving the problem of finding the highest- probability subset under a DPP. For symmetric DPPs, there exists a standard greedy MAP inference algorithm that is linear in M. In this work, we develop version of this algorithm that is also linear for low-rank NDPPs. The low-rank requirement is unique to NDPPs, and highlights the fact that the transformation of the algorithm from the symmetric to the nonsymmetric space is non-trivial. To the best of our knowledge, this is the first MAP algorithm proposed for NDPPs.
We combine the above contributions through experiments that involve learning NDPP kernels and applying MAP inference to these kernels to do subset selection for several real-world datasets. These experiments demonstrate that our algorithms are much more scalable, and that the new kernel decomposition matches the predictive performance of the decomposition from prior work.
2 BACKGROUND
Consider a finite set y= {1, ..) of cardinality M, which we will also denote by [M]. A DPP on [M]] defines a probability distribution over all of its 2 subsets. It is parameterized by matrix L E TR. MMM, called the kernel, such that the probability of each subset Y S [M]] is proportional to the determinant of its corresponding principal submatrix: Pr(Y) x det(Ly, The normalization constant for this distribution can be expressed as a single M x M determinant: Σvgw det(Ly) = (L + T) (Kulesza et al., 2012, Theorem 2)... Hence, Pr(Y)= det(Ly) (onLL))/e det L I). We will use P to denote this distribution.
For intuition about the kernel parameters, notice that the probabilities of singletons { and { are proportional to Li and Ljj. respectively. Hence, is common to think of L's diagonal as representing item qualities. The probability of a pair {i.}} } is proportional to det (t,)) = L L,₃ Thus, if Li < O, this indicates and j interact negatively. Similarly, if L,, L > 0. then i and j interact positively. Therefore, off-diagonal terms determine item interactions. (The vague term "interactions" can be replaced by the more precise term "correlations" if we consider the DPP's marginal kernel instead; see Gartrell et al. (2019, Section 2.1) for extensive discussion.)
In order to ensure that Py defines a probability distribution, all principal minors of L must be non-negative: dy) Z O. Matrices that satisfy this property are called Po-matrices (Fang, 1989, Definition 1). There is no known generative method or matrix decomposition that fully covers the space of all Po matrices, although there are many that partially cover the space (Tsatsomeros, 2004).
One common partial solution is to use a decomposition that covers the space of symmetric Po matrices. By restricting to the space of symmetric matrices, one can exploit the fact that semidefinite (PSD)* (Prussing, 1986). Any symmetric PSD matrix can be written as the Gramian matrix of some set of vectors: L = VVT, where V E R Mx Hence, the VVT decomposition provides an easy means of generating the entire space of symmetric Po matrices. It also has nice intuitive interpretation: we can view the i-th row of V as a length-K K feature vector describing item i.
Po if L is positive
Unfortunately, the symmetry requirement limits the types of correlations that DPP can capture. A symmetric model able to capture only nonpositive interactions between items, since L,ş Lyg =
L3 Z 0. whereas nonsymmetric can also capture positive correlations. (Again, see Gartrell et al. (2019, Section 2.1) for more intuition.) To expand coverage to nonsymmetric matrices in Po. it is natural to consider nonsymmetric PSD matrices. In what follows, we denote by PJ the set of all nonsymmetric (and symmetric) PSD matrices. Any nonsymmetric PSD matrix is in Po (Gartrell et al., 2019, Lemma 1), so Po C Po- However, unlike the symmetric case, the set of nonsymmetric PSD
'Recall that matrix L E MMM is defined he PSD if and only if Læ 2 O, for all e RM
2

---

Published as a conference paper at ICLR 2021
matrices does not fully cover the set of nonsymmetric Po matrices. For example, consider
5/3
= 1/2 1 1 with (1)( det (2)), det( 2 0, but x "Lx < 0 for x
Still, nonsymmetric PSD matrices cover a large enough portion of the Po space to be useful in practice, as evidenced by the experiments of Gartrell et al. (2019). This work covered the Po A, with == VVT for V E RMxK CB for B, C E Mg This decomposition makes use of the fact that any matrix L can be decomposed uniquely as the sum of a symmetrie matrix S = (L + L² and skew-symmetric matrix A = (L )/2. All skew-symmelric matrices A are trivially PSD, since Az = 0 for all a E RM Hence, the L here is guaranteed to be PSD simply because its S uses the standard Gramian decomposition VVT.
space by using the following decomposition: L == and A : BC
In this work we will also only consider Po and leave to future work the problem of finding tractable ways to cover the rest of Po. We propose a new decomposition of that also covers the Po space, but allows for more scalable learning. As in prior work, decomposition has inner dimension K that could be as large as M, but is usually much smaller in practice. Our algorithms work well for modest values of K. In cases where the natural K larger (e.g.. natural language processing), random projections can often be used to significantly reduce K (Gillenwater et al., 2012a).
3 NEW KERNEL DECOMPOSITION AND SCALABLE LEARNING
Prior work on NDPPs proposed a maximum likelihood estimation (MLE) algorithm (Gartrell et al., 2019). Due to that work's 's particular kernel decomposition, this algorithm had complexity cubic in the number of items M. Here, we propose kernel decomposition that reduces this to linear in M.
We begin by showing that our decomposition covers the space of Po matrices. Before diving in,
0 X₂
let us define Σ; = -x 0 as shorthand for 2 x 2 block matrix with zeros on-diagonal and opposite values off-diagonal. Then, our proposed decomposition as follows:
S A, with S VV and A = BCB",
where V. B RRKK and C R is block-diagonal matrix with some diagonal blocks of the form Σ,, with A; > 0, and zeros elsewhere. The following lemma shows that this decomposition covers the space of Po matrices.
Lemma 1. Let A E MMM be skew-symmetric matrix with rank & < M. Then, there exist B E Raxxt xt and positive numbers A
\1//2) A1,...,\(e/21) such that A = BCB where C E Rext is the
block-diagonal matrix with [e/2] diagonal blocks of size 2 given by Σ,, i = 1,..., Le/2] and zero elsewhere.
The proof of Lemma and all subsequent results can be found Appendix F. With this decomposition in hand, we now proceed to show that it can be used for linear-time MLE learning. To do so, we must show that corresponding NDPP ligliikllhood objective and gradient can be computed in time linear in M. Given collection of observed subsets (Yı,.., Yn) composed of items from ) [M], the full formulation of the regularized log- ligeliiildooo is:
(V B. B)C)
log VV + By, CBy. log det (VVT BCB + +1) R(V, B),
where Vy; E RIY|IK denotes matrix composed of the rows of V that correspond to the items in Y, The regularization term, R(V, B), is defined as follows:
M
M
R(V. B) - Σ
+ M
1
(3)
where μ4 counts the number of occurrences of item i in the training set, vs and b, are rows of V and B, respectively, and a; B > 0 are tunable hyperparameters. This regularization is similar to that of prior works (Gartrell et al., 2017; 2019). We omit regularization for C.
3

---

Published as a conference paper at ICLR 2021
Theorem 1 shows that computing the regularized log-likelihood and its gradient both have time complexity linear in M. The complexities also depend on K, the rank of the NDPP, and K', the size of the largest observed subset in the data. For many real-world datasets and we set K = K" Hence, linearity in M means that we can efficiently perform learning for datasets with very large ground sets, which is impossible with the cubic-complexity L decomposition in prior work (Gartrell et al., 2019).
observe that K' < M,
Theorem 1. Given an NDPP with kernel L = VV + BCB", parameterized by V of rank K, B of rank K, and a K K matrix C, we can compute the regularized log- sellikklihoodd (Eq. 2) and its gradient in ( [M + K +nk time, where K' is the size of the largest of the n training subsets.
4 MAP INFERENCE
After learning an NDPP, one can then use it to infer the most probable item subsets in various situations. Several inference algorithms have been well-studied for symmetric DPPs, including sampling (Kulesza & Taskar, 201 I; Anari et al., 2016; Li et al., 2016; Launay et al., 2018; Gillenwater et al., 2019; Poulson, 2019; Derezinski, 2019) and MAP inference (Gillenwater et al., 2012b; Han et al., 2017; Chen et al., 2018; Han & Gillenwater, 2020). We focus on MAP inference:
argmax dey such that =|=k, YCy
(4)
for cardinality budget l: < M. MAP inference is a better fit than sampling when the end application requires the generation of single output set, which is usually the case in practice (e.g., this is usually true for recommender systems). MAP inference for DPPs is known to be NP-hard even in the symmetric case (Ko et al., 1995; Kulesza al., 2012). For symmetric DPPs, one usually approximates the MAP via the standard greedy algorithm for submodular maximization (Nemhauser al., 1978). First. we describe how efficiently implement this for NDPPs. Then, in Section 4.1 we prove a lower bound on its approximation quality. To the best of our knowledge, this is the first investigation of how to apply the greedy algorithm NDPPs.
Greedy begins with an empty set and repeatedly adds the item that maximizes the marginal gain until the chosen set is size k. Here, we design an efficient greedy algorithm for the case where the NDPP kernel is low-rank. For generality, in what follows we write the kernel as L BCB" T, since one can easily rewrite our matrix decomposition (Eq. 1), as well as that of Gartrell et al. (2019), to take this
0
V B
form. For example, for our decomposition: L VV BCB (V B
Using Schur's 's determinant identity, we first observe that, for Y C [M] and E [M], the marginal gain of a NDPP can be written as
dct(Lyy(()) Lis - det(Ly)
= bB By) Cb
(5)
where b, E IRIK and By E R A naïve computation of Eq. 5 is O(K²
2 + k=2), since we must invert YY| x |Y| matrix, where |Y| k. However, one compute Eq. more efficiently by observing that its BB CBC BBy component can actually be expressed without an inverse, as rank-|y | matrix, that be computed in O(K2) time.
Lemma 2. Given B E M CE CeRK and {ai,...,,ay) ax) S [M], let b E RRIk be the i-th row in B and By E RM be matrix containing rows in B indexed by Y. Then, it holds that
BJ BBB Epja
(6)
j=1
where row vectors Pj, qi E x for = 1,...,k satisfy = ba,/(bo, (b,,,CTTT
CbT,). ). qi
bay and
bu, b,, CT 5 q: d.P.
j
Pj+1 = bag C(ba, - ba, CT /
9i+1 = bo, ba, C >
p. qi-
(7)
Σ T²
i=1
4

---

Published as a conference paper at ICLR 2021
Algorithm I Greedy MAP inference/conditionigg for low-rank NDPPs
I: Input: B I RM C E R x the cardinality
D And {aI,.... 0k for conditioning
2: initialize + []. Q f
and Y +
3: A; + -b, for i E [M] where b, E RIIK is the i-th row in B
4: & argmax; A, and Y (a)
aa t aj for conditioning
5: while |Y] < k do
6:
p (ba - by P C)
7:
q+ - b, - bbCCP
8:
P P and Q + [Q:q]
9:
A A, (b,CpT) (b,C q ) for i [M]igy
argmax, A and Y aa
f- pat
a
d|Y|+1 for conditioning
11: end while
refurn Y
D return {A, for conditioning
Table 1: Algorithm complexities for several DPP models. Our model and the symmetric DPP model (Gartrell al., 2017) can perform both tasks time linear the size ground set M, but ours is more general model that can capture positive as well as negative item correlations.
MLE Leaming
MAP Inference
MLE Learning
MAP Inference
Low-rank DPP Models
Runtime
Runtime
Memory
Memory
Symmetric DPP (Gartrell et al., 2017)
O(MK + a
O + MK?)
O(MK)
O(MK)
Nonsymmetric DPP (Gartrell et al., 2019)
O(M + MK² nK*)
O(MKK MK*)
O(M2)
O(MKY'
Scalable nonsymmetric DPP (this work)
O(MK² nK*)
O(MKK MK3
O(MK K²)
O(MK K
Plugging Eq. 6 into Eq. 5, the marginal gain with respect to Y U {a} can be computed by simply updating from the previous gain with respect to Y. That is,
det det (Lrutao))))) b. det det((r)((a))
|++1
Σ (b.Cr)) (b,CT,))) F=1
(8)
del(Lrutt) - (D.CPTVH)) (b.C g det(Lr)
The marginal gains when Y = 0 are equal to diagonals of L and require O(MK²) operations. Then, computing the update terms in Eq. for all i E [M]] needs O(MK) operations. Since the total number of updates is h. the overall complexity becomes O(MK2 + MKk). We provide a full description of the implied greedy algorithm for low-rank NDPPs in Algorithm
Table summarizes the complexitiy of our methods and those of previous work. Note that the full M x + matrix used compute the DPP normalization constant in Gartrell et al. (2019), which is why this approach has memory complexity of O(M2) for MLE learning.
4.1 APPROXIMATION GUARANTEE FOR GREEDY NDPP MAP INFERENCE
As mentioned above, Algorithm an instantiation of the standard greedy algorithm used for submodular maximization (Nemhauser al., 1978). This algorithm has (1 app)-opproiiiiaooo guarantee for the problem of maximizing nonnegative, monotone submodular functions. While the function f(Y) log det (Ly) is submodular for a symmetrie PSD L (Kelmans Kimelfeld, 1983), it is not monotone. Often, as in Han & Gillenwater (2020), it is assumed that the smallest eigenvalue of L greater than which guarantees montonicity. There is no particular evidence that this assumption is true for practical models, but nevertheless the greedy algorithm tends to
The exact memory complexity for MAP inference is 3M since V, B, and C used in this model are all
M K matrices.

---

Published as a conference paper at ICLR 2021
perform well in practice for symmetric DPPs. Here, we prove similar approximation guarantee that covers NDPPs as well, even though the function f(Y) = log det(Ly) is non-submodular when L is nonsymmetrie. In Section 5.5, we further observe that, as for symmetric DPPs, the greedy algorithm seems to work well in practice for NDPPs.
We leverage recent result of Bian et al. (2017), who proposed extension of greedy algorithm guarantees to non-submodular functions. Their result based on the submodularity ratio and curvature of the objective function, which measure to what extent it has submodular properties. Theorem 2 extends this to provide an approximation ratio for greedy MAP inference of NDPPs.
Theorem 2. Consider a nonsymmetric low-rank DPP VV +
BCB", where ,,B B.B are of rank K, and C & RR* x Given cardinality budget k. let Juin and max denote the smallest and largest singular values of Ly for all Y S [M] and |Y] < 2k. Assume that min > 1. Then,
- e(1---1/4)
log det (dee(LLy)) 2 2(log O maxx /log miin) log det ())
(10)
where YG is the output of Algorithm and Y the optimal solution of MAP inference Eq. 4.
Thus, when the kernel has small value of log maxx log mmin, the greedy algorithm finds near- optimal solution. In practice, we observe that the greedy algorithm finds a near-optimal solution even for large values of this ratio (see Section 5.5). As remarked above, there is no evidence that the condition min > 1 is usually true in practice. While this condition can be achieved by multiplying L by constant, this leads to a (potentially large) additive term in Eq. 10. We provide Corollary in Appendix D, which excludes the o min
assumption, and quantifies this additive term.
4.2 GREEDY CONDITIONING FOR NEXT-ITEM PREDICTION
We briefly describe here a small modification to the greedy algorithm that is necessary if one wants to use it tool for next-item prediction. Given a set Y [M] Kulesza al. (2012) showed that a DPP with L conditioned on the inclusion of the items in Y forms another DPP with kernel L := Ly-LLyLL Ly.y Ly- Ly.y where Y [M]VY: The singleton probability Pr(Yu(i} Y) L can be useful for doing next-item prediction. We can use the same machinery from the greedy algorithm's 's marginal gain computations to effectively compute these singletons. More concretely, suppose that we are doing next-item prediction as a shopper adds items to a digital cart. We predict the item that maximizes the marginal gain, conditioned on the current cart contents (the set Y). When the shopper adds the next item to the cart, we update Y to include this item, rather than our predicted item (line 10 in Algorithm 1). We then iterate until the shopper checks out. The comments on the righthand side of Algorithm 1 summarize this procedure. The runtime of this prediction is the same that of the greedy algorithm, O(MK² + (MKY)) We note that this cost is comparable to that of an approach based on the DPP dual kernel from prior work (Mariet et al., 2019), which has OMM + K3 + |)] complexity. However, since is non-trivial to define the dual kernel for NDPPs, the greedy algorithm may be the simpler choice for next-item prediction for NDPPs.
EXPERIMENTS
To further simplify learning and MAP inference, we set B V, which results in L VV VCV = V + CVV This change also simplifies regularization, so that we only perform regularization on V, as indicated in the first term of Eq. leaving us with the single regularization hyperparameter of a. While setting B = V restricts the class of nonsymmetric L kernels that can be represented, we compensate for this restriction by relaxing the block diikoniaaooall structure imposed on C, so that we learn a full skew-symmetric K x K matrix C. To ensure that C' and thus A is skew-symmetric, we parameterize C' by setting C = D D". were D varies over K
Code for all experiments is available at https:///e http:///ccoooooggtteeeeeeeeeeeee able- www.-ooommmmtiiiiieeeee
5.1 DATASETS
We perform experiments several real-world public datasets composed of subsets:

---

Published as a conference paper at ICLR 2021
1. Amazon Baby Registries: This dataset consists of registries "baskets" of baby products, and has been used in prior work on DPP learning (Gartrell et al., 2016; 2019; Gillenwater et al., 2014; Mariet & Sra, 2015). The registries contain items from 15 different categories, such as "apparel", with a catalog of up to 100 items per category. Our evaluation mirrors that of Gartrell et al. (2019); evaluate the popular apparel category, which contains 14,970 registries, well as on dataset composed of the three most popular categories: apparel, diaper, and feeding, which contains a total 31,218 registries.
2. UK Retail: This dataset (Chen et al., 2012) contains baskets representing transactions from an online retail company that sells unique all-occasion gifts. We omit baskets with more than 100 items, leaving us with a dataset containing 19,762 baskets drawn from a catalog of M = 3,941 products. Baskets containing more than 100 items are in the long tail of the basket-size distribution of the data, so omitting larger baskets is reasonable, and allows us to use a low-rank factorization of the DPP with K 100.
3. Instacart: This dataset (Instacart, 2017) contains baskets purchased by Instacart users. We omit baskets with more than 100 items, resulting in 3.2 million baskets and a catalog of 49,677 products.
4. Million Song: This dataset (McFee et al., 2012) contains playlists ("baskets"" of songs played by Echo Nest users. We trim playlists with more than 150 items, leaving with 968,674 baskets and catalog of 371, 371,410 songs.
5.2 EXPERIMENTAL. SETUP AND METRICS
We use small held-out validation set, consisting of 300 mandmlly-selctted baskets, for tracking convergence during training and for tuning hyperparameters. A random selection of 2000 of the remaining baskets are used for testing, and the rest are used for training. Convergence is reached during training when the relative change in validation log- oogllikeihhood is below predetermined threshold. We use PyTorch with Adam (Kingma & Ba, 2015) for optimization. We initialize C from the standard Gaussian distribution with mean 0 and variance I, and B (which we set equal to V) is initialized from the uniform(0, 1) distribution.
Subset expansion task. We use greedy conditioning to do next-item prediction, described in Section 4.2. We compare methods using a standard recommender system metric: mean percentile rank (MPR) (Hu al., 2008; Li et al., 2010). MPR of 50 is equivalent random selection; MPR of 100 means that the model perfectly predicts the next item. See Appendix A for a complete description of the MPR metric.
Subset discrimination task. We also test the ability of a model to discriminate observed subsets from randomly generated ones. For each subset in the test set, generate a subset of the same length by drawing items uniformly at random (and we ensure that the same item is not drawn more than once for a subset). We compute the AUC for the model on these observed and random subsets, where the score for each subset is the log- likeelhooodd that the model assigns to the subset.
5.3 PREDICTIVE PERFORMANCE RESULTS FOR LEARNING
Since the focus of our work is on improving NDPP scalability, we use the low-rank symmetric DPP (Gartrell et al., 2017) and the low-rank NDPP of prior work (Gartrell et al., 2019) baselines for our experiments. Table 2 compares these approaches and our scalable low-rank NDPP. We see that NDPPs generally outperform symmetric DPPs. Furthermore, we that scalable NDPP matches or exceeds the predictive quality of the baseline NDPP. We believe that our model sometimes improves upon this baseline NDPP due to the use of simpler kernel decomposition with fewer parameters, likely leading to simplified optimization landscape.
5.4 TIME COMPARISON FOR LEARNING
In Fig. I, we report the wall-clock training time the decomposition of Gartrell et al. (2019) (NDPP) and our scalable NDPP for the Amazon: 3-category (Fig. I(a)) and UK Retail (Fig. 1(b)) datasets. As expected, we observe that the scalable NDPP trains far faster than the NDPP for datasets with large ground sets. For the Amazon: 3- 3-category dataset, both approaches show comparable results, with the scalable NDPP converging 1.07 x faster than NDPP. But for the UK Retail dataset, which has much larger ground set, scalable NDPP achieves convergence about 8.31 x faster. Notice
7

---

Published as a conference paper at ICLR 2021
Table 2: Average MPR, AUC, and test log-likelihood for all datasets, for the low-rank symmetric DPP (Gartrell et al., 2017), low-rank NDPP (Gartrell et al., 2019), and our scalable NDPP models. MPR and AUC results show 95% confidence estimates obtained via bootstrapping. Bold values indicate improvement over the symmetric low-rank DPP outside of the confidence interval. See Appendix B for the hyperparameter settings used in these experiments. The baseline NDPP model cannot be feasibly trained on the Instacart and Million Song datasets, as memory and computational costs are prohibitive due to large ground set sizes.
Amazon: Apparel (M = 100)
Amazon: 3-category (M = 300)
Metric
Sym 62.63 + 1811
Nonsym
Scalable nonsym
Sym
Nonsym
Scalable nonsym
MPR
72.20 + 3.07
69.02 2.57
610 = 2.73
74.10 2.49
73.04 2.58
AUC
0.68 = 0.05
0.77 = 0.03
0.74 + 0.03
0.76 = 0.03
0.82 + 0.02
0.82 = 0.02
test log-likelhhood
-10.02
-9.64
-9.63
-18.11
-1699
-17.14
UK Retail (M = 3,911)
Instacart (M = 49,677)
Metrie
Sym
Nonsym
Scalable nonsym
Sym
Nonsym Scalable nonsym
MPR
69.95 1.32
74.17 4 1.37
76.79 + 1.17
93.86
0.55
93.13 F 0.53
AUC
0.58 + 0.01
0.66 + 0.01
0.73 0.01
0.83 0.01
0.85 0.005
test log- likelikoeiood
-116.23
104.38
-100.65
-72.81
-72.74
Million Song (M
371,410)
Metric
Sym
Nonsym
Scalable nonsym
MPR AUC
90.37 = 0.71
90.41 : 0.75 0.77 = 0.01
0.69 0.01
test likelikdl
-335.25
-317.16
negative log-liken likelihood
200
negative log- ligliiel
600
calabbe NDPP NDPP
Scalable
NDPP
150
400
100
200
50
200
400
1000
2000
wall clock time (sec)
wall clock time (sec)
(a) Amazon: 3-category
(b) UK Retail
Figure The negative log-likelihood of the training set for Gartrell et al. (2019)'s NDPP (blue, dashed) and our scalable NDPP (red, solid) versus wall-clock time for the (a) Amazon: 3-category and (b) UK Retail datasets.
that our scalable NDPP also opens to the door to training on datasets with large M, such as the Instacart and Million Song dataset, which is infeasible for the baseline NDPP due to high memory and compute costs. For example, NDPP learning using Gartrell et al. (2019) for the Million Song dataset would require approximately 1.1 TB of memory, while using our scalable NDPP approach requires approximately 445.9 MB.
5.5 PERFORMANCE RESULTS FOR MAP INFERENCE
We run various approximatation algorithms for MAP inference, including the greedy algorithm (Algo- rithm 1), stochastic greedy algorithm (Mirzasoleiman al., 2015), MCMC-based DPP sampling (Li et al., 2016), and greedy local search (Kathuria & Deshpande, 2016). The stochastic greedy algorithm computes marginal gains of a few items chosen uniformly at random and selects the best among them. The MCMC sampling begins with a random subset Y of size k and picks i E Y and j & Y uniformly at random. Then, it swaps them with probability det det(Lrug)een and iterates this process. The greedy local search algorithm (Kathuria & Deshpande, 2016) starts from the output from the greedy algorithm, Y and replaces E y with j & y that gives the maximum improvement, if such i,j exist. This replacement process iterates until improvement exists, or
+ det(Ly)) (Ly))
8

---

Published as a conference paper at ICLR 2021
Table 3: Average relative enror and 95% confidence intervals of MAP inference algorithms on NDPPs learned from real-world datasets. For all datasets, we evaluate 10 kernels learned with different initializations, and run 100 random trials for stochastic greedy (Mirzasoleiman al., 2015) and MCMC sampling (Li et al., 2016). All emror are relative to greedy local search (Kathuria & Deshpande, 2016).
Algorithms
Amazon: Apparel
Amuzon: 3-cutegory
Instacurt
Million Song
Greedy (Algorithm I)
0.0336 + 0.0666
0.0093 + 0.0015
0.0446 + 0.M05
0.0173 = 0.0028
0.0052 = 0.0017
Stochustic greedy
0.1606 + 0.0133
0.1838 + 0.0116
0.0960 0.0078
0.1229 = 0.0091
0.0823 0.0108
MCMC sampling
0.7155 + 0.0287
0.7094 + 0.0207
0.9365 + 0.0342
1.9291 0.047
1.0493 0.0607
Table 4: Wall-clock time (in milliseconds) of MAP inference algorithms on NDPPs learned from real-world datasets.
Algorithms
Amazon: Apparel
Amazon: 3-category
UK Retuil
Instacart
Million Song
Greedy local search
5.78 ms
9.67
58.74
1.024
7.277
Greedy (Algorithm
ms
0.34 ms
1.60 ms
36.16 ms
338.09 ms
Stochustic greedy
0.25 ms
0,47
1.79 ms
36.94 ms
348.67 ms
MCMC sampling
0.35
2.85 ms
42.85 ms
303.20 ms
most h- (10k) steps have been completed, to guarantee a tight approximation (Kathuria & Desh-
pande, 2016). We use greedy local search as a baseline since it always returns a better solution than greedy. However, it is the slowest among all algorithms, its time complexity is O(M KK* log k). We choose & = 10, and provide more details of all algorithms in Appendix C.
To evaluate the performance of MAP inference, report the relative log-determinant ratio defined as
log (oettLy)--oe det Ly)
log det(Lv.)
where Y is the output of benchmark algorithms and Y* is the greedy local search result. Results are reported in Table 3. We observe that the greedy (Algorithm 1) achieves performance close to that of the significantly more expensive greedy local search algorithm, with relative errors of up to 0.045. Stochastic greedy and MCMC sampling have significantly larger errors.
For completeness, in Appendix E we also present experiments comparing the performance of greedy and exact MAP on small synthetic NDPPs, for which the exact MAP can feasibly computed.
5.6 TIME COMPARISON FOR MAP INFERENCE
We provide the wall-clock time of the above algorithms for real-world datasets in Table 4. Observe that the greedy algorithm is the fastest method for all datasets except Million Song. For Million Song, MCMC sampling faster than other approaches, but it has much larger relative errors in terms of log-determinant (see Table 3), which is not suitable for our purposes.
6 CONCLUSION
We have presented a new decomposition for nonsymmetric DPP kernels that can be learned in time linear in the size of the ground set, which is significant improvement over the complexity of prior work. Empirical results indicate that this decomposition matches the predictive performance of the prior decomposition. We have also derived the first MAP algorithm for nonsymmetric DPPs and proved lower bound on the quality of its approximation. In future work we hope to develop intuition about the meaning of the parameters in the matrix and consider kernel decompositions that cover other parts of the nonsymmetric Po space.
9

---

Published as a conference paper at ICLR 2021
REFERENCES
Nima Anari, Shayan Oveis Gharan, and Alireza Rezaci. Monte Carlo Markov Chain Algorithms for Sampling Strongly Rayleigh Distributions and Determinantal Point Processes. In Conference on Learning Theory (COLT), 2016.
Andrew An Bian, Joachim M. Buhmann, Andreas Krause, and Sebastian Tschiatschek. Guarantees for Greedy Maximization of Non-submodular Functions with Applications. In International Conference on Machine Learning (ICML), 2017.
Erdem Bıyık, Kenneth Wang. Nima Anari, and Dorsa Sadigh. Batch Active Learning Using Determi- nantal Point Processes. arXiv: 1900:099067775 2019.
Victor-Emmanue Brunel. Learning Signed Determinantal Point Processes through the Principal Minor Assignment Problem. In Neural Information Processing Systems (NeurIPS), 2018.
Daqing Chen, Sai Laing Sain, and Kun Guo. Data mining for the online retail industry: case study of RFM model- based customer segmentation using data mining. Journal of Database Marketing & Customer Strategy Management, 2012.
Laming Chen, Guoxin Zhang, and Eric Zhou. Fast greedy MAP inference for Determinantal Point Process to improve recommendation diversity. Neural Information Processing Systems (NeurlPS), 2018.
Michat Dereziński. Fast determinantal point processes via distortion-fres intermediate sampling. In Conference on Learning Theory (COLT), 2019.
Mohamed Elfeki, Camille Couprie, Morgane Riviere, and Mohamed Elhoseiny. GDPP: Learning Diverse Generations using Determinantal Point Processes. In International Conference on Machine Learning (ICML), 2019.
Li Fang. On the Spectra of P- and Po-Matrices. In Linear Algebra and its Applications, 1989.
Mike Gartrell, Ulrich Paquet, and Noam Koenigstein. Bayesian low-rank determinantal point processes. In Conference on Recommender Systems (RecSys). ACM, 2016.
Mike Gartrell, Ulrich Paquet, and Noam Koenigstein. Low-Rank Factorization of Determinantal Point Processes. In Conference on Artificial Intelligence (AAAI), 2017.
Mike Gartrell, Vietor-Emmanuel Brunel, Elvis Dohmatob, and Syrine Krichene. Learning Non- symmetric Determinantal Point Processes. In Neural Information Processing Systems (NeurIPS),
2019.
Jennifer Gillenwater, Alex Kulesza, and Ben Taskar. Discovering Diverse and Salient Threads in Document Collections. In Empirical Methods in Natural Language Processing (EMNLP), 2012a.
Jennifer Gillenwater, Alex Kulesza, and Ben Taskar. Near-Optimal MAP Inference for Determinantal Point Processes. In Neural Information Processing Systems (NIPS), 2012b.
Jennifer Gillenwater, Alex Kulesza, Emily Fox, and Ben Taskar. Expectation- ExpectationMaximizattie for learning Determinantal Point Processes. In Neural Information Processing Systems (NIPS), 2014.
Jennifer Gillenwater, Alex Kulesza, Zelda Mariet, and Sergei Vassilvtiskii. A Tree-Based Method for Fast Repeated Sampling of Determinantal Point Processes. In International Conference on Machine Learning (ICML), 2019.
Insu Han and Jennifer Gillenwater. MAP Inference for Customized Determinantal Point Processes via Maximum Inner Product Search. In Conference Artificial Intelligence and Statistics (AISTATS), 2020.
Insu Han, Prabhanjan Kambadur, Kyoungsoo Park. and Jinwoo Shin. Faster greedy MAP inference for determinantal point processes. In International Conference on Machine Learning (ICML), 2017.
10

---

Published as a conference paper at ICLR 2021
Yifan Hu, Yehuda Koren, and Chris Volinsky. Collaborative Filtering for Implicit Feedback Datasets.
In International Conference on Data Mining (ICDM), 2008.
Instacart. The Instacart Online Grocery Shopping Dataset, 2017. URL https: setps://wwet /www.
instacart. com/ /ttp///aeeeeeeeeeeeeeoarreeeeeee
2017. Accessed May 2020.
Tarun Kathuria and Amit Deshpande. On sampling and greedy map inference of constrained determinantal point processes. arXiv preprint 6001v:160..01551,222 2016.
A.K. Kelmans and B.N. Kimelfeld. Multiplicative submodularity of a matrix's principal minor as function of the set of its rows and some combinatorial applications. Discrete Mathematics, 1983.
Diederik P Kingma and Jimmy Ba. Adam: A method for stochastic optimization. In International Conference on Learning Representations (ICLR), 2015.
Chun-Walie Ko, Jon Lee, and Maurice Queyranne. An Exact Algorithm for Maximum Entropy Sampling. Operations Research, 1995.
Alex Kulesza and Ben Taskar. Learning determinantal point processes. In Conference on Uncertainty in Artificial Intelligence (UAT), 2011. 1.
Alex Kulesza, Ben Taskar, et al. Determinantal Point Processes for Machine Learning. Foundations and Trends@ in Machine Learning, 2012.
Claire Launay, Bruno Galerne, and Agnès Desolneux. Exact Sampling of Determinantal Point Processes without Eigendecomposition. arXiv preprint arXiv:/802.0844: 2018.
Chengtao Li, Stefanie Jegelka, and Suvrit Sra. Fast DPP Sampling for Nystrom with Application to Kernel Methods. In International Conference on Machine Learning (ICML), 2016.
Yanen Li, Jia Hu, ChengXiang Zhai, and Ye Chen. Improving One-class Collaborative Filtering by Incorporating Rich User Information. In Conference on Information and Knowledge Management (CIKM), 2010.
Zelda Mariet and Suvrit Sra. Fixed-point algorithms for learning Determinantal Point Processes. In International Conference on Machine Learning (ICML), 2015.
Zelda Mariet and Suvrit Sra. Diversity Networks: Neural Network Compression Using Determinantal Point Processes. In International Conference on Learning Representations (ICLR), 2016.
Zelda Mariet, Mike Gartrell, and Suvrit Sra. Learning Determinantal Point Processes by Sampling Inferred Negatives. In Conference on Artificial Intelligence and Statistics (AISTATS), 2019.
Brian McFee, Thierry Mahtiu- Daniel PW Ellis, and Gert RG Lanckriet. The million song dataset challenge. In Proceedings of the 21 International Conference World Wide Web, 2012.
Baharan Mirzasoleiman, Ashwinkumar Badanidiyuru, Amin Karbasi, Jan Vondrák, and Andreas Krause. Lazier Than Lazy Greedy. In Conference on Artificial Intelligence (AAAT), 2015.
G. Nemhauser, L. Wolsey, and M. Fisher. An Analysis of Approximations for Maximizing Submodu- lar Set Functions I. Mathematical Programming, 14(1), 1978.
Jack Poulson. High-performance sampling of generic Determinantal Point Processes. arXiv preprint arXiv:1905.00166 00165, 2019.
John E. Prussing. The Principal Minor Test for Semidefinite Matrices, Journal of Guidance, Control, and Dynamics, 1986.
Aidean Sharghi, Ali Borji, Chengtao Li, Tianbao Yang, and Boqing Gong. Improving Sequential Determinantal Point Processes for Supervised Video Summarization. In Proceedings of the European Conference on Computer Vision (ECCV), 2018.
G Thompson. Normal forms for skew symmetycmt matrices and Hamiltonian systems with first integrals linear in momenta. In Proceedings of the American Mathematical Society, 1988.
11

---

Published as a conference paper at ICLR 2021
Robert C Thompson. Principal submatrices IX: Interlacing inequalities for singular values of submatrices. Linear Algebra and its Applications, 1972.
Michael Tsatsomeros. Generating and Detecting Matrices with Positive Principal Minors, In Focus on Computational Neurobiology, 2004.
Mark Wilhelm, Ajith Ramanathan, Alexander Bonomo, Sagar Jain, Ed H. Chi, and Jennifer Gillen- water. Practical Diversified Recommendations on YouTuue with Determinantal Point Processes. In Conference on Information and Knowledge Management (CIKM), 2018.
Cheng Zhang, Hedvig Kjellström, and Stephan Mandt. Determinantal Point Processes for Mini-Bate Batch Diversification. In Conference Uncertainty in Artificial Intelligence (UAI), 2017.
12

---

Published as a conference paper at ICLR 2021
A MEAN PERCENTILE RANK
We begin our definition of MPR by defining percentile rank (PR). First, given a set J, let Pt.J UPJJUO) J). The percentile rank of an item given a J is defined
PR..J Ever 1(p,J 2 pe.J) 100%
DJJ
where V indicates those elements in the ground set y that are not found in J.
For our evaluation, given test set Y, select random element i and compute y We then average over the set of all test instances T to compute the mean percentile rank (MPR):
MPR
Σ YET
PR.Y
T
B HYPERPARAMETERS FOR EXPERIMENTS TABLE 2
Preventing numerical instabilities: The first term on the right side of Eq. 2 will be singular whenever [Yl > K, where Y, is an observed subset. Therefore, to address this in practice we set to the size of the largest subset observed in the data, K', in Gartrell al. (2017). However, this does not entirely address the issue, as the first term on the right side of Eq. 2 may still be singular even when |Y|| S K. In this case though, we know that we are not at maximum, since the value of the objective function is -c. Numerically, to prevent such singularities, in our implementation we add small el correction each Ly, when optimizing Eq. (we set -5 in our experiments).
We perform grid search using a held-out validation set to select the best performing hyperparameters for each model and dataset. The hyperparameter settings used for each model and dataset are described below.
Symmetric low-rank DPP (Gartrell al., 2016). For this model, we use K for the number of item feature dimensions for the symmetric component V, and for the regularization hyperparameter for V. We use the following hyperparameter settings:
Both Amazon datasets: K = 30, 0.
UK Retail dataset: = 100, 1.
Instacart dataset: K = 100, a = 0.001.
Million Song dataset: K = 150, a 0.0001.
Baseline NDPP (Gartrell et al., 2019). For this model, to ensure consistency with the notation used Gartrell al. (2019), we use 1) denote the number of item feature dimensions for the symmetric component V, and D' to denote the number of item feature dimensions for the nonsymmetric components, B and C. As described in Gartrell et al. (2019), a the regularization hyperparameter for the V. while 3 and are the regularization hyperparameters for B and C, respectively, We use the following hyperparameter settings:
Both Amazon datasets: D 30, 0.
Amazon apparel dataset: D' 30.
Amazon three-category dataset: D' = 100.
UK Retail dataset: D = 100, D' = 20, a 1.
All datasets: 3 =
Scalable NDPP. As described in Section 3, use K to denote the number of item feature dimensions for the symmetric component V and the dimensionality the nonsymmetric component C. a the regularization hyperparameter. We use the following hyperparameter settings:
Amazon apparel dataset: K = 30, a = 0.
Amazon three-category dataset: K 100, 1.
UK dataset: K 100, 0.01.
13

---

Published as a conference paper at ICLR 2021
Instacart dataset: K = 100, a 0.001.
Million Song dataset: K = 150, 0.01.
For all of the above model configurations we use a batch size of 200 during training, except for the scalable NDPPs trained on the Amazon apparel, Amazon three-category, Instacart, and Million Song datasets, where a batch size of 800 is used.
C BENCHMARK ALGORITTMS FOR MAP INFERENCE
We test the following approximate algorithms for MAP inference:
Greedy local search. This algorithm starts from the output of greedy, Y and replaces i y with j É y that gives the maximum improvement of the determinant, if such i, exist. Kathuria & Deshpande (2016) showed that running the search for such a swap O(k2log(H log(k/e)) times with an accuracy parameter gives tight approximation guarantee for MAP inference for symmetric DPPs. We set the number of swaps to [A_ log(10k)] for E = 0.1 and use greedy local search as baseline, since is strictly an improvement on the greedy solution. The proposed greedy conditioning can be used for fast greedy local search. Specifically, for each E y Algorithm I can compute marginal improvements conditioned by yc \ {i} in time OMk,, and thus its runtime can be O(M KIAAA log(k/e)). However, it is the slowest among all of our benchmark algorithms.
Stochastic greedy. This algorithm computes marginal gains of a few items chosen uniformly at random and selects the best among them. Mirzasoleiman et al. (2015) proved that (M/k) log(1/e) samples are enough to guarantee an (1- 1/e a apappoooinat ratio for submodular functions (i.e., symmetric DPPs). We choose this setting, the time complexity of stochastic greedy is (MMK" log(1/E)). which is better than the naïve exact greedy algorithm. However, implement (A (Algorithm 1). This is because the stochastic greedy uses different random samples for every iteration and this does not take advantage of the amortized computations in Lemma 2. In our experiments, we simply modify line 10 in Algorithm for stochastic greedy (argmax is operated on random subset of marginal gains), hence can run in O(M (M/k) log(1/s)) time. In practice, we observe that stochastic greedy is slightly slower than exact greedy due to the additional costs of the random sampling process.
0.1 and set the number of samples to (M/k) log(10)]. Under
note that it is worse than that of efficient greedy
MCMC sampling. We also compare inference algorithms with sampling from a nonsymmetric DPP. To the best of our knowledge, exact sampling of a non-Hermitian DPP was studied in Poulson (2019), which requires the Cholesky decomposition with ((M)) ) complexity. This infeasible for large M. To resolve this, Markov Chain Monte-Carde (MCMC) based sampling is preferred (Li et al., 2016) for symmetric DPPs. In particular, we consider a Gibbs sampling for k-DPP, which begins with random subset Y with size k. and picks e Y and j & Y uniformly at random. Then, it swaps them with probability
det (et(Lru())))
(1I)
det(Lyupryr)) + dy)
and repeat this process for several steps. Li et al. (2016) showed that O(Nk log(k/E)) swaps are enough to approximate the ground-truth distribution under symmetric DPPs. However, for a fair runtime comparison to Algorithm 1, we set the number of swaps to [3N/K].
D COROLLARY OF THEOREM 2
Theorem requires the technical condition aain > but in practice there is particular evidence that this condition holds. While this condition can be achieved by multiplying L by a constant, this leads to a (potentially large) additive term in Eq. 10. Here, we provide Corollary which excludes the Jmin > assumption from Theorem 2, and quantifies this additive term.
Corollary 1. Consider a nonsymmetric low-rank DPP L= VV BCB where V, B are of rank K, and C E R.K x Given cardinality budget k. let a and max denote the smallest and
14

---

Published as a conference paper at ICLR 2021
approxination ratio
approximation ratio
singular values
10f
0.4
0.4
0.2
0.2
-... Thoorem
---
1/44 guarandoon
10
10
10
10
20
log moox/ mmin)
log maxx/ mmin)
index
(a) Symmetric DPP
(b) Nonsymmetric DPP
(c) Singular values
Figure 2: Approximation ratios of greedy with respect to different values of log(O m from Corollary under (a) symmetric DPP and (b) nonsymmetric DPP. (c) The singular values of the kernels learned for the "Amazon: 3-category" dataset. We construct 10,000 random Po matrices L E TR with rank K = 3, whose singular values are from the learned kernels.
largest singular values of Ly for all C [M] and [У] < 2k. Let K. == Tmaxx F Then,
e 1/4)
4(1 4(1---1/41
log det degdleLaa))
2 K + log det (olet(Ly.)
1o (1 log miin) (12)
2 2 Bogk+* +1
where yG the output of Algorithm and Y the optimal solution of MAP inference Eq. 4.
The proof of Corollary is provided in Appendix F.5. Note that instead of Corollary has a log(O max/ min)) term in the denominator.
log((riin)
E PERFORMANCE GUARANTEE FOR GREEDY MAP INFERENCE
The matrices learned on real datasets are too large to compute the exact MAP solution, but we can compute exact MAP for small matrices. In this section, explore the performance of the greedy algorithm studied in Theorem for 5 X synthetic kernel matrices. More formally, we first pick K singular values $1, $2, $3 from kernel learned for the "Amazon: 3-category" dataset (a plot of these singular values can be seen in Fig. 2(c)) and generate L = Vidiag([s1, $2, S3 where Vi, V/2 E Rx are random orthonormal matrices. To ensure that L is a Po matrix, we repeatedly sample Vi, V22 until all principal minors of L are nonnegative. We also evaluate the performance of the symmetric DPP, where the kernel matrices are generated similarly to the NDPP, except we set We set k = 3 and generate 10,000 random kernels for both symmetric DPPs and NDPPs.
Vi V2.
The results for symmetric and nonsymmetric DPPs are shown in Fig. 2(a) and Fig. 2(b), respectively, We plot the approximation ratio of Algorithm 1, i.e., log (Ly )/ log det Ly with respect to m,n)), from Corollary 1. We observe that the greedy algorithm for both often shows approximation ratios close to 1. However, the worst-case ratio for NDPPs is worse than that of symmetric DPPs; log det(Ly for L E Pt is non-submodular, and the greedy algorithm with a nonsubmodular function does not have tight of worst-case bound as in the symmetric case.
15

---

Published as a conference paper at ICLR 2021
F PROOFS
F.1 PROOF OF LEMMA
Lemma 1. Let A E RMM be skew-symmetric matrix with rank & < M. Then, there exist B R&xt and positive numbers (a...,,/E/22) such that A BCB where C Rxxt is the block-diagonal matrix with [e/2] diagonal blocks of size given by E,, I,..., Le/2] and zero elsewhere.
Proof. First, we note that rank of a nonsingular skew-symmetric matrix is always even, because all of its eigenvalues are purely imaginary and come in conjugate pairs. There exists some orthogonal matrix P E RMM xM and
0
X1
A1
0
0
X2
A2
0
Σ
(13)
0
Xe/2
0
XE/2
0
0
such that = P (see, (Thompsnn,,,n, 1988, Proposition 2.1)).
1:
Let C be the x supmatrix of obtained by keeping its first rows and columns and let Q =
0
where I is the l x é identity matrix. Then, Σ = QCQ and one can write A = PQCQ P Setting B PQ proves the lemma.
F.2 PROOF OF THEOREM
Theorem 1. Given NDPP with kernel L VV + BCB", parameterized by V of rank K, B of rank K, and a K x K matrix C, we can compute the regularized log-likee log-likelihood (Eq, 2) and its gradient in OMK?2 K +k33)) time, where K' is the size of the largest of the training subsets.
Proof. We first show that the log-likelihood can be computed in time linear in M. Using the matrix determinant lemma, one can easily verify that the DPP normalization term can be computed as
det(I + L) det (I+ BC)
)
V B
det [2x +
(V
BC)
(14)
B
where Ia is the identity matrix with dimension 2K. As Eq. 14 requires matrix- catrii--mutiplliationn between (2K) x M matrices and the determinant of (2K) x (2K) matrices, this allows us to transform O(13) operation into O(MK² + K3) one.
Having established that the normalization term the likelihood can be computed O(M + K3) time, we proceed with characterizing the complexity of the other terms in the likelihood. The first term in Eq. consists of determinants of size [Y]- Assuming that these never exceed size K', each be computed in most ((K*)) time. The regularization tenm is simple sum of norms that can be computed O (M)) time. Therefore, the full regularized log-likelihond can computed in O(M R2 + K3 33 time.
16

---

Published as a conference paper at ICLR 2021
To prove that the gradient of the log- legliike can be computed in time linear in M, begin by showing that the logarithm of DPP normalization term can be factorized as follows:
Z
log dot(I I + L)
(15)
Ik ))
log det
12K
(V
B)
(16)
B
0
(or
IK
= log det
(V
B)
+ log det
(17)
C-1
B
IK
VTV
V
log det
+ log dct(C)
(18)
V
C-I+Bl
= log det
x + VTV) + log det (C- + B (I V(Ix + +VV + log det (C)
(19)
where Eq. 17 follows from the determinant commutativity (i.e., det(AB) = det (A) det(B)) and Eq. 18 and Eq. 19 come from the Schur's determinant identity For simplicity, write
X = -- I + VV- and (C-1) = =C-T, and note that X depends only on V.
The gradient of Z has three parts: VZ ( V VZZ) where each can be computed as
VvZ = Vv log det I + +V + Vv logdettt B
(20)
+ V)-
-XB((C- B XB) (C- + BB) B BV XV
(21)
BB VB log det XB)
(22)
= -B ((C-1 + XB) + +C- +B XB)
(23)
OZ = Vc leg Vc log (C+ B XB)
(24)
= = = -C- -C-
XB) C-T
(25)
Observe that X combines a M x M identity matrix with M x K matrices, hence multiplying it with a M K matrix (e.g., XV or X XB) can be computed O(M)2) time. Since each of the remaining matrix inverses in Eq. 21, Eq. 23, and Eq. 25 involve K x K matrix inverse, with cost of ((K3) operations, we have net computational cost O(MK² + K²) for computing V log del (I + L).
The gradient of the first term in Eq. 2 involves computing gradients of determinants of size at most 2 (log (et(()) = (A ")
K', which results in size K' matrix inverses, since for matrix A, JA.,
Each of these inverses can be computed in O(K*3) time. The gradient of the simple sum-of-norms regularization term can be computed in O(MK) time. Therefore, combining these results with the results above for the complexity of the gradient of the normalization term, we have the following overall complexity of the gradient for the full log-likelihood: O(M + K3 + nK*3).
F.3 PROOF OF LEMMA 2
Lemma 2. Given B E RMMK C E RXKK, and = [aj...., ax} S. [M], let b, E RRIk be the
i-th row in B and By E | be a matrix containing rows in B indexed by Y. Then, it holds that
k Σr}
(By -By
p3 li
(6)
1=1
where row vectors Pj, q, RRIX for 1,...,k satisfy Pı /ar/ CbT). di = =bay and
bu, ba, CT 5 q.P:
gi+1 bog ba; C Σ
Pj+1
p: gi.
bos C(ba, bas CT )/ 9: Tp)
ddt C D
= (et(A)det((D) -CA
17

---

Published as a conference paper at ICLR 2021
Proof. We prove by induction on k. When &o the result is trivial because
(By By b,, (ba, Cb ( pi Q1-
(26)
Now we assume that the statement holds for k Let Y = {a1.-.., ; 1) and a == dg. From the inductive hypothesis, it holds
ax--))
h-1
(By -By ΣpJar
(27)
/=1
Now write
Byutat (B -1
Byula)
(28)
= -Brula)
((xx)
by
C CB}
b7))
Bru
(29)
ByyBBy ba BB
By baC
-1
By
= -By
b
(30)
ba
To handle the inverse matrix we employ the Schur complement, which yields
X
y
0)
X X-1yz y-=y-X-- -x-ty
+ w - - -x-1y y
-x
(31)
for any non-singular square matrix X E Ø column vector y E R and row vector E I2Ixk unless (w- - 2X -1y) = 0. Applying this, we have
By
Bb 1 bB bb = (B,CB)) 9)
b. - bB By -
(By (B,CB b, CB (B,CB ) Cb
bb (By
1.
(32)
Substituting Eq. 32 into Eq. 30, we obtain
Byuta (Bruta) CB CByuta) Bru(a}
(33)
B (BrCBT) -1 By + (b - CBB (b, (,
b.C By BB,
(34)
k-1
(b k- p. CCC) (b. - 13
Σ Pj +
(35)
baC b. bT Σ P,
1=1
k-1
Σ p. 1=1
+ pk qs
(36)
where the third line holds from the inductive hypothesis Eq. 27 and the last line holds from the definition of Pk q R lxx
F.4 PROOF OF THEOREM 2
Theorem 2. Consider nonsymmetric low-rank DPP L - VV + BCB where V B are of rank K, and C E x Given a cardinality budget lel amin and Fmax denote the smallest and largest singular values of Ly for all Y S [M] and |Y| < 2k. Assume that a >> 1. Then,
4(1 -114)
log det (deetLLya))))) 2 2(log maax maxx//og mmin)
log (ct((L-))) (Ly.)
(10)
where G is the output of Algorithm and Y* is the optimal solution of MAP inference in Eq. 4.
18

---

Published as a conference paper at ICLR 2021
Proof. The proof of Theorem relies on approximation guarantee for nonsubmodular greedy maximization (Bian et al., 2017, Theorem 1). We introduce their result below.
Theorem 3 ((Bian al., 2017, Theorem 1)). Consider set function f defined all subsets of (l....................... ..] = [M]] is monotone nondecreasing and nonnegative, i.e., 0 S f(Y) S f(X) for VY X X C [M]. Given a cardinality budget k 2 1, let Y* be the optimal solution of ma f(Y) and y 0, y = {al,...,0e))
1.....k. k be the successive chosen by the greedy algorithm with budget k. Denote yy be the largest scalar such that
(UY' U (()) - f(X")) 2U((XUY)-
(37)
for _ [M]] |X| and O,...,k - and be the smallest scalar such that
S(Y*-1 U - S(Y'-I(U)))) 2 (1 a) u f(Y-))))
(38)
for C [M]. X = k and i E @ X. Then, it holds that
s(y*) 1 (1-- f(Y*).
(39)
In order to apply this result for MAP inference of NDPPs, the objective should be monotone nondecreasing and nonnegative. We first show that min 1 is sufficient condition for both
monotonicity and nonnegativity.
Lemma 3. Given Po matrix L E: RMM and the budget k 0, senftctioon f(Y) log dyy defined on Y S [M] is monotone nondecreasing and nonnegative for [Y] S when m >
The proof of Lemma 3 is provided in Appendix F.6. Next, we aim to find proper bounds and y. To resolve this, we provide the following upper and lower bounds of the marginal gain for f(Y) log det d.
Lemma 4. Let f(Y) = log det (Ly) and assume that amin > 1. Then, for Y S [M], Y < 2k and i Y, it holds that
f(Y U (1) f(Y) 2 log miin-
(40)
J(Y U (i)) - J(Y) < 2 log Tmax -
llg Fmin
(41)
where min and max are the smallest and largest singular values of Ly for all Y S [[]]] |Y| < 2k.
The proof of Lemma is provided in Appendix F.7. To bound y, we consider X C [M]] X and denote X Y .er + 0. Then,
Σ (UYY - ((Y)) Σ ( - S(Y*) 2 llg Tmin
(42)
VEXIV
where the last inequality comes from Eq. 40. Similarly, we get
(xur))-(YY)) =
u
(43)
1=1
I <r(2 log maax llog mmin)
(44)
where the last inequality comes from Eq. 41. Combining Eq. 42 to Eq. 44, we obtain that
Etexyy. f(X* f(X*U {it) - f(X*)
log Fmin
(45)
f(XUY9- - -(Y't
2 log mmax log amin
which allows to choose 7
log - 2
Iog Gwin
To bound a, we similarly use Lemma to obtain
u - ((XUY))) f(y-) f(Y
log Cinin
(46)
log Tmax log o
19

---

Published as a conference paper at ICLR 2021
log 2loog Omax - log Omina) log Omox log Emin
and we can choose log Omax
log anio
Now let
log Piuis log Canin
Then 7
and
- 2k
Putting y and into Eq. 39, we have
1 - -1 a
2 2x -1
(
2(w-1) (2n-1)2
1 -
(47)
2(K 1)
2x
e - 2(K - 1)
2 2(K 1)
(48)
(2x 1)2
4 exxp( 1/4)
(49)
2x 1
where the second inequality holds from the fact that max,21 (2x-1)2
2(x-1)
4 and 1 e_
4 4exp(-1/444 for æ E [0, 1/4].
F.5 PROOF OF COROLLARY
Corollary 1. Consider a nonsymmetric low-rank DPP L - BCB where V are of rank K, and C E BKKK Given cardinality budget k, let min and mmxx denote the smallest and largest singular values of Ly for all S [M] and [Y] S 2k. Let K. = Tmax/ T Then,
4(1 - -1/4)
1(1 (11---1/1 1-11))
log det (eddet(Lyy Lya) 2
log (eed(e((-)))
2log + 21og*+1
ke (1 - log minn) (12)
2log +1
where yo the output of Algorithm and Y the optimal solution of MAP inference Eq. 4.
Now consider L' -- where e is the exponential constant. Then, Juin Jmin Ominn e and d'nax = Fmax (- Using the fact that log det(L(y) = log dee(Lyy
Proof.
|Y] log Jmin, we obtain the result.
F.6 PROOF OF LEMMA 3
Before stating the proof, we introduce interlacing properties of singular values.
Theorem (Interlacing Inequality for Singular Values, (Thompson, 1972, Theorem 1)). Consider A
real matrix A E Mx with singular values 0ı 2 02 2 ... Min((M and its supmatrix
B T PxQ with singular values B1 2 B2 2 2 Pmin(P.Q) Then, the singular values have the following interlacing properties:
0 2 3r,
for 1, .,... min(P.))
(PQQ):
(50)
3: Z
for i == l.... min( P P M.Q-MMP+QL N).
(51)
Note that when M N and P Q N 1, it holds that Br 2 0or = 1,..., N - 2.
We are now ready to prove Lemma
Lemma 3. Given a Po matrix L E MxM and the budget 2 0, set function f(Y) = log det Ly) defined on Y C [M] is monotone nondecreasing and nonnegative for [Y] S k when o > 1.
Proof. Since L Po. all of its principal submatrices also in Po. By the definition of a Po matrix, it holds that
[de)(Lr)) = det (e = II o. (Lr)
(52)
where (Lv) for i E [[Y]] are singular values of Ly. nonnegative for all Y such that |Y| I K due to ox (Ly) amin 1. Similarly, we have
Then, F(Y) = E, log(a, (Ly)) is
F(Y U {a}) - F(Y) = I1 log o, (Lruça)) - EM log o, (Ly) 2 log J > 0 from Eq. 50.
20

---

Published as a conference paper at ICLR 2021
F.7 PROOF OF LEMMA 4
Lemma 4. Let S(Y) = log det (v) and assume that mmin > 1. Then, for Y C [M], Y| < 2k and
i Y, it holds that
( 2 log min:
(40)
f(Y U (t)) - f()) < 2 log maax - log miin
(41)
where min and max are the smallest and largest singular values of Ly for all Y [M]] |Y|
2k.
Proof. For Po matrix, remark that its determinant is equivalent to the product of all singular values. For Y S [M] and & Y, from the interlacing inequality of Eq. 50 we have that
I+++ Σ
F
log o'
Σ log og 2 log 0 log mmin f=1
(53)
where of and oy are the j-th largest singular value of Lru(i) and Ly, respectively. Similarly, using Eq. 51, we get
F(Y U (+)) F(Y) log(o(o') -log011 S 2 log mux
log Onin-
(54)
21


Reviews:
[{"review_id": "mLQ656HfVPt", "reviewer": "AnonReviewer2", "summary": "", "questions": "", "limitations": "", "rating": 9, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "This paper presented a novel kernel decomposition for nonsymmetric determinantal point processes, which enables linear time of inference and learning w.r.t. the cardinality of the ground set M. This is a significant improvement over previous arts and makes NDPP practical in relatively large datasets. The theoretical of the paper is solid and supportive to the main claim of the paper. This paper is well written and easy to follow. Even for readers without much theoretical background, this paper is still moderately friendly since the logic and insight are clear. I vote this paper a strong acceptance, except for some minor concerns as follows:\n\n1) The authors employed a way to simplify the kernel decomposition below Theorem 1. However, it is unclear what the impact of this simplification is to the exactness of learning or inference. It would me a plus if the authors can give some theoretical analysis on the gap between such simplification and $P_0^+$.\n\n2) I notice the authors provide both learning and inference procedures for NDPP. It seems these two procedures can formulate a way to learn latent NDPP in an Expectation-Maximization fashion. I would suggest the authors give some discussion about the feasibility of such integration, and this can be an interesting direction.\n\nIn general, I enjoy reading this paper and think this paper is insightful.", "extra_scores": {}, "extra_text": {}}, {"review_id": "KIWKldKHUA", "reviewer": "AnonReviewer4", "summary": "", "questions": "", "limitations": "", "rating": 7, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "Nonsymmetric determinantal point processes (NDPPs) received some attention recently because they allow modeling of both negative and positive correlations between items. This paper developed scalable learning and MAP inference algorithms with space and time complexity linear in ground set size, which is a huge improvement compared to previous approaches. Experimental results show that the algorithms scale significantly better, and can roughly match the predictive performance of prior work.\n\nThis is a well written paper and I recommend its acceptance. Scalable learning and MAP inference algorithms are important for the application of the NDPPs model, which seems promising compared with its symmetric counterpart in experiments.\n\nI have some (minor) comments listed below.\n\n1. In Lemma 1, the result is only proved for skew-symmetric matrices with even rank. Does it hold for odd rank matrices? This is important to support the claim that the new decomposition covers the $P_0^+$ space.\n\n2. Equation (3) uses notation $\\lambda_i$, which is already used in Lemma 1. This could cause confusion.\n\n3. In the paragraph after Theorem 1, it is proposed to set $B$ = $V$ and relax $C$. Is this used in Section 4? If not, I would suggest moving it to the experiments section, and adding some comparison in Table 2 to show the impact of this simplification.\n\n4. I cannot quite understand the last sentence before Lemma 2. What can be computed in $O(K^2)$ time?\n\n5. The footnote in Table 1 might cause confusion because it can be mis-interpreted as a square.\n\n6. In G.1, the first sentence after equation (13), do you mean when $M$ is odd or when $\\ell$ is odd?\n\n7. In equation (24), $X$ should be $B^TXB$\n\n8. The inverse of $C$ appears in the gradient of $Z$. Is $C$ guaranteed to be invertible in the learning algorithm? And how are $V$, $B$, $C$ initialized in the algorithm?\n\n9. In equation (31), please double check if we need the reciprocal in the denominator.\n", "extra_scores": {}, "extra_text": {}}, {"review_id": "54b_th9wRo", "reviewer": "AnonReviewer3", "summary": "", "questions": "", "limitations": "", "rating": 8, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "This paper propose a decomposition for non-symmetric determinantal point process (NDPP) kernels (M*M) which reduces the requirements of storage and running to linear in cardinality (M). Additionally, they derive a NDPP maximum a posteriori inference algorithm that applies to both their proposed kernel and the previous work (NDPP). In their experiments, they show both learning kernels and the MAP inference for subset selection on real-life datasets. \n\nPros:\n\n\t○ This paper is well-written and easy to follow.\n\t○ The author provide sufficient calculation process and relevant proofs of scalable NDPP.\n\t○ For the existing problems of traditional DPP method is time consuming, need large memory and could not apply to large set, the scalable NDPP really solves them (e.g., this method could run on Instacart and Million Song datasets). I think it is practical. \n\t\n\nMajor Concern:\n\n\t○ I have some doubts about the authenticity of the experimental results in Table 2, for the following reason: in the previous work, i.e., [1] Learning Nonsymmetric Determinantal Point Processes. Gartrell et al. Neurips2019. The results under average MPR have signifcant difference with same datasets and same hyperparameter settings. However, in NDPP paper,  for Amazon: Apparel, the MPR of sym dpp is 77.42±1.12, MPR of nonsym DPP is 80.32±0.75, for uk retail,  the MPR of sym dpp is 76.79±0.6, MPR of nonsym DPP is 79.45±0.57.  In this paper,   for Amazon: Apparel, the MPR of sym dpp is 62.63±1.81, MPR of nonsym DPP is 72.2±3.07, for uk retail,  the MPR of sym dpp is 69.95±1.32, MPR of nonsym DPP is 74.17±1.37. There is a huge gap btw these two versions. I would like to know what causes this gap.", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: This paper proposes a technique of decomposing the nonsymmetric kernel of determinantal point processes, which enables inference and learning in time and space linear with respect to the size of the ground set.  This substantially improves upon existing work.  The proposed method is well supported both with theory and experiments.  All of the reviewers find that the contributions are significant, and no major flaws are identified through reviews and discussion.  The determinantal point process might not be one of the most popular topics in the ICLR community today but certainly is relevant.
Decision: Accept (Oral)