Learning Multinomial Logits in $O(n\log n)$ Time

A multinomial logit (MNL) model supported on the universe $U=[n]=\{1,\mathrm{\dots },n\}$ is specified by a set $\{w_1,\mathrm{\dots },w_n\}$ of $n$ positive values called weights. A slate is a non-empty subset of $[n]$. An MNL $M$, for any given slate $S\subseteq [n]$, induces a conditional distribution denoted $M_S$ whose support is $S$ and where the probability of each item $i\in S$ is given by:¹¹1The terminology we adopt comes from the Economics literature [42], this is called a logit model because if we let $w_i=e^{{\theta }_i}$, then $M_S(i)=\text{softmax}{(S)}_i=e^{{\theta }_i}/{\sum }_{j\in S}{e}^{{\theta }_j}$.

An MNL $M$ can be accessed by a $\text{Sample}$ oracle, which operates as follows: given a slate $S$, $\text{Sample}(S)$ returns $i\in S$ chosen according to the distribution $M_S$. Given MNLs $M$ and $M^{\prime }$, we define two notions of distance between them:

In this paper we obtain algorithms that approximate an unknown MNL $M$ in $d_1$ distance, while our lower bounds apply even to the less challenging problem of obtaining estimates with small $d_{\mathrm{\infty }}$ distance.

In this paper, we study algorithms for the MNL learning problem. We obtain two algorithms, one using adaptive queries and the other using non-adaptive queries. Both algorithms query only slates of size two and run in time proportional to their query complexity. Our adaptive algorithm makes $O\left(\frac{n}{{\epsilon }^3}\log n\right)$ queries; we give a lower bound of $\mathrm{\Omega }(\frac{n}{{\epsilon }^2}\log n)$ queries. Summarizing:

Our non-adaptive algorithm makes $O\left(\frac{n^2}{{\epsilon }^3}\log n\log \frac{n}{\epsilon }\right)$ queries; this is complemented by a lower bound of $\mathrm{\Omega }(\frac{n^2}{{\epsilon }^2}\log n)$. Summarizing:

As we mentioned above, the lower bounds described also hold for the weaker $d_{\mathrm{\infty }}$ distance.

Our results are surprising: for a constant $\epsilon$, our seemingly harder problem can be solved as fast as (noisy) sorting. Furthermore, our lower bounds hold for unit-time oracle queries of any slate size. Hence, restricting the algorithms to slates of size two incurs no loss in efficiency.

Existing methods, especially ones that approximate the winning distribution over the universal slate $[n]$, do not seem to apply to our problem. As a simple example of the difficulty, consider an algorithm that guarantees an ${\mathrm{\ell }}_1$-estimate of the full slate distribution within an error of $\epsilon \in (0,1/2)$. Consider now the MNL on $\{1,2,3\}$ with weights $w_1=1-\epsilon$, $w_2=w_3=\epsilon /2$. Suppose the algorithm returns the estimate ${\widehat{w}}_1=1-\epsilon$, ${\widehat{w}}_2=\frac{3\epsilon }{4}$, ${\widehat{w}}_3=\frac{\epsilon }{4}$; clearly, ${\parallel w-\widehat{w}\parallel }_1=\frac{\epsilon }{2}\le \epsilon$. But, $\left|\frac{w_2}{w_2+w_3}-\frac{{\widehat{w}}_2}{{\widehat{w}}_2+{\widehat{w}}_3}\right|\ge \frac{1}{4}$, and therefore the algorithm cannot guarantee small error on the slate $\{2,3\}$. Similarly, prior work that additively estimates the winning distributions on all size-two slates cannot be used to obtain a good approximation on every slate. Moreover, noisy sorting algorithms also do not imply our guarantees. Indeed, while our final algorithm will perform a noisy sorting step at the beginning, the key challenge lays in performing the additional comparisons that allow the accurate estimate of each conditional distribution, as we detail below.

It is not difficult to obtain a cubic time algorithm for our problem. Indeed, consider the naive algorithm that works as follows. For each pair $\{i,j\}$ of items in the universe, estimate $w_i/w_j$ to within a $(1\pm \epsilon )$-multiplicative error, or declare that their ratio (or its inverse) is larger than $\frac{n}{\epsilon }$. One can easily show that this algorithm will need to query each pair $\approx \frac{n\log n}{{\epsilon }^3}$ times to guarantee these bounds; the total cost would then be $\approx \frac{n^3\log n}{{\epsilon }^3}$. From the output of this algorithm, it is easy to approximate the output distribution for any slate.³³3Indeed, if the items of this slate have weights that are within a $\frac{n}{\epsilon }$ factor of each other, the $(1\pm \epsilon )$-approximation error will make it possible to approximate the winning probability of any item to within a $(1\pm O(\epsilon ))$-factor (so that the total variation error will be at most $O(\epsilon )$). If, instead, the slate contains pairs $\{i,j\}$ of items such that , then the lighter item $i$ will have a probability of winning in the slate not larger than $O(\frac{\epsilon }{n})$, and hence we can estimate its winning probability to be zero – given that there are at most $n-1$ such light items in a slate, the total variation error is no larger than $O(\epsilon )$.

With some effort, this algorithm can be improved. The idea is to carefully control the pairs of items, querying only pairs that are nearby in the order induced by the weights (which can be approximately inferred via a noisy sorting algorithm). To avoid querying too many nearby pairs, one can first cluster items whose weights are within a constant factor, in a query-efficient manner, and then select a center from each cluster. One can determine the weight ratio of each item to its cluster center, and the weight ratios of successive items in the sorted list of cluster centers. This method can be shown to produce an algorithm that compares $O(n)$ pairs of items, with each such comparison performing $\approx \frac{n}{{\epsilon }^3}{\log }^{3}n$ queries. While this yields a quadratic time algorithm, it is unclear how this can be further improved to being quasi-linear.

We construct our quasi-linear adaptive algorithm by building on the clustering idea described above. We first partition the universe into clusters of similar-weight items and select a center for each cluster. We then estimate the ratio of the weights $w_i/w_c$ for every item $i$ in the cluster with center $c$. Finally, we construct a forest on the cluster centers, where every edge is labeled with an estimate of the ratio between the weights of the two centers it connects. We call this data structure the estimation-forest. This allows us to obtain estimates of the ratio of the weights for arbitrary pairs of items by combining these estimates along paths in the forest.

Following this strategy, the error compounds multiplicatively along the paths. To circumvent this issue without requiring more accurate ratio estimates – which would lead to a higher complexity – we design the forest so that any two centers whose weights the algorithm might want to compare are at a short distance from each other. To achieve this property, the topology of the forest is constructed adaptively.

Additionally, to further improve the sample complexity (and runtime), we dynamically adjust the number of queries required to approximate the weight ratio between two centers. In particular, if the total weight of items lighter than a given item $i$ is not large enough with respect to the weight of $i$, then it becomes unimportant to estimate the ratio of the weight of $i$ over the weight of any of these items.

Our estimation-forest data structure dynamically determines query sequences for weight-ratio estimation and enables all cluster-center–cluster-center and cluster-item–cluster-center comparisons in $O\left(\frac{n}{{\epsilon }^3}\log n\right)$ queries.

While the above algorithm is adaptive, we also obtain a non-adaptive version. The idea is to first design a new adaptive algorithm that queries each pair of items only $O\left(\frac{1}{{\epsilon }^3}\log n\log \frac{n}{\epsilon }\right)$ times. We then query every pair a fixed number of times and then simulate this new adaptive algorithm on the precomputed answers; this leads to a non-adaptive algorithm with $O\left(\frac{n^2}{{\epsilon }^3}\log n\log \frac{n}{\epsilon }\right)$ queries.

The two lower bounds in our paper are proved using a reduction from the problem of identifying the $k$ coins with the highest heads probability in a collection of $n$ biased coins.

In particular, we construct an MNL supported on an even-sized universe whose items are divided into pairs, each pair representing the two sides of a coin. We then order the pairs so that the weights of the items in a pair are much larger than those in preceding pairs. This way, we can assume without loss of generality that any learning algorithm is only querying slates corresponding to our original pairs.

In Section 2 we review related work. Section 3 introduces key tools and notation that we use throughout the paper. Section 4 gives an overview of our adaptive algorithm, while Section 5 gives an overview of our non-adaptive algorithm. Section 6 gives an overview of our adaptive and non-adaptive lower bounds. Finally, in Section 7 and Section 8, we conclude with some open questions.

All the proofs and algorithms missing from the paper can be found in the full version.

A large body of the literature focuses on finding MNL weights maximizing the likelihood of a collected dataset. In this setting, usually, the queries are either fixed [45, 19, 15, 31] or sampled from a distribution [32, 29, 30, 27]. When the dataset actually comes from a hidden MNL model, some of these algorithms guarantee that the estimated (normalized) weights approximate the hidden (normalized) weights [29, 30, 27, 37, 39, 36].

These works are not directly applicable to our setting because of the following two main issues.

1. (i)

   Approximately recovering the normalized weights is equivalent to providing a good estimate of the winning distribution of the full slate. However, this is insufficient to accurately estimate the winning distribution for smaller slates, as we discussed in the Introduction.

2. (ii)

   Most of these works assume that the maximum ratio between two weights is upper bounded by a constant [29]. Note that such an assumption would greatly simplify our problem given that we could accurately estimate the weights with respect to any anchor item. Therefore, the interesting setting is one where there is no a priori bound on the ratio of the weights. Furthermore, as these algorithms are non-adaptive, they are subject to our lower bound of $\mathrm{\Omega }(n^2\log n)$ queries for our problem.

Stepping outside the task of fitting the weights themselves, Falahatgar et al. [17] adaptively query $O\left(\frac{n\cdot \log (n)\cdot \min \{n,1/\epsilon \}}{{\epsilon }^2}\right)$ slates of size two (i.e., pairs) and produce an additive estimate within $\epsilon$ for all other pairs, in a family of models that are more powerful than MNLs [44]. However, in order for an additive approximation of the slates of size two to generalize to all other slates with an ${\mathrm{\ell }}_{\mathrm{\infty }}$-error of ${\epsilon }^{\prime }$, one needs $\epsilon \le \frac{{\epsilon }^{\prime }}{n}$. Therefore, applying their algorithm as a blackbox would require $\mathrm{\Omega }(n^4\log n)$ queries. Moreover, their algorithm heavily relies on providing an additive approximation – specifically, each estimated probability is rounded to the closest multiple of $\epsilon$. Hence, it appears hard to generalize their work to larger slates even in a non-blackbox manner.

Other works considered minimizing the ${\mathrm{\ell }}_{\mathrm{\infty }}$-difference of the logarithm of the weights, or approximating the weight ratios of all the pairs [43, 20, 40, 24]. These guarantees would be sufficient to obtain a good estimate for all slates, but are not necessary. Indeed, if $w_1=e^{1000}$, $w_2=e^{\epsilon }$ and ${\widehat{w}}_1=e^{100}$, ${\widehat{w}}_2=e^{\epsilon }$, then ${\parallel \ln (w)-\ln (\widehat{w})\parallel }_{\mathrm{\infty }}$ can be made arbitrarily large by increasing the constant $1000$, yet, $\widehat{w}$ provides good predictions: $|\frac{w_1}{w_1+w_2}-\frac{{\widehat{w}}_1}{{\widehat{w}}_1+{\widehat{w}}_2}|\le \epsilon$. Moreover, these works analyze non-adaptive algorithms and assume a bound on the ratio of the weights, and therefore suffer our $\mathrm{\Omega }(n^2\log n)$ lower bound.

We also mention a separate line of work that focuses on developing statistical tests to determine whether a given dataset of comparisons is consistent with an MNL model [38, 25, 34]. These works are complementary to ours in that they address model validation rather than estimation, and they do not provide algorithms for learning the underlying MNL parameters.

Other classical problems involving MNLs include: (i) sorting/ranking the weights [18, 17, 9, 41, 35], (ii) finding the top-$k$ items with largest weight [21, 10, 11, 12, 22], and (iii) finding the item of maximum weight [16, 26]. Some works make assumptions about the weights (e.g., adjacent weights are sufficiently separated) and seek an exact output with high probability [21], while others do not make further assumptions but only require a probably approximately correct (PAC) output [41]. These problems have also been explored in the “dueling bandits” literature [3].

While we will use an $O(\frac{n\log n}{{\epsilon }^2})$ approximate sorting algorithm by Falahatgar et al. [17] as a first step in our algorithm, these results are not sufficient by themselves to learn an MNL under our definition. Our lower bound will be proved by showing that the top-$\frac{n}{2}$ problem (and the ranking problem) reduces to our MNL learning problem. Interestingly, despite this, we obtain an $O(\frac{n\log n}{{\epsilon }^3})$ learning algorithm that is only $O(\frac{1}{\epsilon })$-factor worse than the best possible algorithm for MNL ranking [17].

Our problem can also be described in the context of conditional sampling. Let $\mu$ be a hidden distribution over $[n]$. Algorithms can, adaptively, make the following types of queries: chosen a set $S\subseteq [n]$, an oracle returns an item of $S$ sampled according to distribution $\mu$ conditioned on $S$.⁴⁴4If $S$ has probability 0, the oracle returns a uniform at random item from $S$. The goal in distribution testing is usually to make the smallest number of queries to establish whether $\mu$ satisfies certain properties, such as, e.g., uniformity [5]. However, the problem of estimating the probability $\mu (i)$, for $i\in [n]$, has also been considered [8, 6, 2].

Our problem, on the other hand, asks for the minimum number of queries to accurately estimate $\mu (i\mid S)$ for each $i\in S\subseteq [n]$.⁵⁵5In our proofs, we will assume that the weights are strictly positive for simplicity. However, the same algorithms also work when weights of zero are allowed. Indeed, the distribution $\mu$ can be seen as the weights of an MNL $M$ and therefore, $\mu (i\mid S)=M_S(i)$ for $i\in S\subseteq [n]$. Observe that having an estimate only for $\mu (i)$ is equivalent to an estimate of the winning distribution on the full slate – insufficient to estimate the winning distribution for smaller slates. Some algorithms provide a multiplicative $(1\pm \epsilon )$ estimate for $\mu (i)$ for $i\notin B$ where $B$ is a set such that ${\sum }_{b\in B}\mu (b)\le \epsilon$; this is a stronger property than an additive approximation of the full slate. However, it still cannot provide accurate estimates for slates that are either subsets of $B$ or that span across $B$ and $[n]\setminus B$. Note also that it can be $|B|=\mathrm{\Theta }(n)$ meaning that the distribution of most slates cannot be estimated. From a technical standpoint, Chakraborty et al. [8] achieve this guarantee by building a complete binary tree where edges are labeled with probabilities. This idea bears some high level similarities with our estimation-forest, but details differ. Indeed, their tree is static while the topology of our forest is adaptively chosen, which is crucial for a tight $O(n\log n)$ bound. Moreover, their tree is populated by querying slates of arbitrary size, while we only query slates of size two. Finally, as argued above, the guarantees provided by their tree are insufficient to estimate the winning distributions of all slates. Recent work has also focused on different query models [1, 33, 28]; however, their results are incomparable to ours.

We remark that a common paradigm for designing distribution testing algorithms in the traditional (unconditional) setting is that of testing by learning (see, e.g. [7]), in which a property is tested by first approximately learning the underlying distribution, and then checking whether the learned distribution has the property in question. Our work serves as a conditional counterpart to this paradigm that works for the more challenging case in which the property being tested requires approximating the behavior of all conditional distributions.

MNLs are a special case of RUMs; hence, algorithms for learning RUMs on all slates could be used to learn an MNL. However, the best known algorithms for general RUM learning require exponentially many queries to slates of size $\mathrm{\Theta }(\sqrt{n})$ [13]. In contrast, we show that MNLs can be learned using only $O(n\log n)$ queries to slates of size two.

An ordered clustering of $[n]$ is given by an ordered partition $(C_1,\mathrm{\dots },C_T)$ of $[n]$ and a corresponding list $(c_1,\mathrm{\dots },c_T)$ of centers such that $c_i\in C_i$ for each $i$. In our algorithms, the centers can always be thought to be sorted in increasing order of weight. Here, for $v\in [n]$, let $\gamma (v)\in [T]$ be the unique index such that $v\in C_{\gamma (v)}$; we call $\gamma (v)$ the cluster index of $v$.

Let $F=([n],E)$ be an undirected forest supported on $[n]$. For $u,v\in [n]$ in the same connected component of $F$, let $P(u,v)$ be the (unique) path in $F$ from $u$ to $v$. We use $d(u,v)$ to denote the (unweighted/hop) distance in $F$ between vertices $u$ and $v$, where if $u$ and $v$ are in different connected components, we define $d(u,v)=\mathrm{\infty }$.

Let $\overrightarrow{E}:=\{(u,v)\in V^2\mid \{u,v\}\in E\}$. A directed weighting of the edges of $F$ is a function $r:\overrightarrow{E}\to {ℝ}_{>0}$ such that $r(u,v)=1/r(v,u)$. For a path $P=u_1,\mathrm{\dots },u_t$ in $F$ define $r(P)={\prod }_{i=1}^{t-1}r(u_i,u_{i+1})$, and if $t=1$, let $r(P)=1$.

Interestingly, for the purpose of constructing $\widehat{M}$, it turns out that the specific value of $t$ is irrelevant. This observation allows us to reduce the problem of learning $M$ to that of constructing a $(t,\epsilon )$-estimation-forest for a single, arbitrary choice of $t$. The central challenge then lies in designing an efficient topology for the estimation-forest.

The most natural topology would be a path on the items (after a noisy-sorting step). However, along a path, two items of comparable weight can be separated by a super-constant distance $d=\omega (1)$. To preserve property 2 of Definition 4, one would then need to construct a $(d,\epsilon )$-estimation-forest, which would incur a query cost of $\mathrm{\Omega }(n{d}^2\log n)$. Since $d$ can be as large as $\mathrm{\Theta }(n)$, this is clearly suboptimal, suggesting the need for a topology with low diameter. Note that even a complete binary tree also can yield super-constant length paths, implying an $\omega (n\log n)$ query cost.

On the other hand, to achieve a very small diameter, one might consider a star or a tree topology with unbounded arity. However, in these cases, one would need to estimate extremely large weight ratios, leading to high query complexity. This, in fact, explains why a disconnected graph is required.

Another natural direction would be to consider general (non-acyclic) graphs. In fact, we could consider a path with skips to decrease the diameter (perhaps exploiting modern shortcutting results [23]). The difficulty is that, in a cyclic graph, the weight of an item depends on the particular path chosen, and different paths can yield inconsistent estimates. Thus, acyclicity of the topology is essential to ensure that an explicit MNL can be extracted from it.

In the full version, we provide an efficient algorithm for constructing an $(O(1),\epsilon )$-estimation-forest. The resulting topology takes the form of a forest of lobster graphs: items are first clustered together, as described below, and a forest of unbounded-arity trees is then constructed over the resulting cluster centers. The arity of each tree is not predetermined but is instead adaptively chosen as the algorithm progresses, in order to balance estimation accuracy and query efficiency. Interestingly, the diameter of the trees in our forest can be super-constant. However, each tree will have the property that if two items are at distance more than $O(1)$, then one of the two is so much larger than the other that their ratio can be taken to be infinite without incurring a large error. Thanks to this property, from the perspective of any single item, one can consider the tree to have constant diameter and lose at most $O(\epsilon )$ in the final estimate.

We now go more into the details of our solution to efficiently build an estimation-forest. When constructing the estimation-forest, some ratio estimates might be costlier to obtain than others. In order to maintain a low query complexity, we leverage the fact that if two items have similar weights, fewer queries are required to estimate the ratio of their weights. In the first step to build our estimation-forest, we exploit this observation via a pre-processing step, which sorts the items in approximately increasing order of weights, and produces clusters of similar items resulting in a cluster graph, defined as follows.

In the full version we provide an algorithm to efficiently build a cluster graph. Our algorithms employs a noisy sorting procedure of Falahatgar et al. [17] as a subroutine and builds on it to partition the vertices and compute the edge weights $r$. We show in Figure 1 a cluster graph produced by our algorithm.

Observe that a cluster graph is not yet an estimation-forest. Indeed, there might be items in different clusters (but close in the ordering) whose ratio is constant. To obtain an estimation-forest, we add extra edges between some pairs of centers. We do so in an iterative way, starting from the center of the last cluster and moving backwards. A priori, these multiplicative estimates can potentially be costly to obtain, since the ratio between the weights of distinct cluster centers could be arbitrarily large. In order to maintain a low query complexity, we employ a careful thresholding strategy. This ensures that we only require an accurate estimate of the ratio when this is not too large to make a significant difference in the MNL winning distributions. For instance, if the ratio between two items is greater than $\mathrm{\Omega }(\frac{n}{\epsilon })$, then it is safe to act as if the second item’s weight is infinitely larger than the first, as this approximation only causes a $d_1$-error of magnitude $O(\epsilon )$. When we find two clusters that are incomparable, we restart the iteration process from the last cluster that was comparable. It can be shown that this leads to an $(O(1),\epsilon )$-estimation-forest. We show in Figure 2 a forest that can be produced by our algorithm.

In summary, our algorithm consists of three phases. In the first phase, we construct a $(\mathrm{\Theta }(1),\mathrm{\Theta }(1),\mathrm{\Theta }(\epsilon ))$-cluster graph. In the second phase, we extend the cluster graph to a $(\mathrm{\Theta }(1),\mathrm{\Theta }(\epsilon ))$-estimation-forest. Finally, in the third phase, we use the forest to recover an estimate of the MNL weights. A representation of the steps in our algorithm is in Figure 3. The first two phases require at most $O(\frac{n\log n}{{\epsilon }^3})$ queries, while the last one does not make any further queries, yielding our main result: