A tractable discrete fractional programming: application to constrained assortment optimization

Abstract

In this note, we consider a discrete fractional programming in light of a decision problem where limited number of indivisible resources are allocated to several heterogeneous projects to maximize the ratio of total profit to total cost. For each project, both profit and cost are solely determined by the amount of resources allocated to it. Although the problem can be reformulated as a linear program with $O(m^2 n)$ variables and $O(m^2 n^2)$ constraints, we further show that it can be efficiently solved by induction in $O(m^3 n^2 \log mn)$ time. In application, this method leads to an $O(m^3 n^2 \log mn)$ algorithm for assortment optimization problem under nested logit model with cardinality constraints (Feldman and Topaloglu, Oper Res 63: 812–822, 2015).

Appendix

1.1 A Generating the individual candidate sets

Here we propose a simple algorithm to generate all $\mathcal {T}_i(c)$ for nest i and $c = 1, \dots , n$. Also, we provide an analysis that, by our generating procedure, $\sum _{C_i = 1}^n \vert \mathcal {T}_i(C_i) \vert = O(n^2)$. This result further tightens the $O(n^3)$ bound in Feldman and Topaloglu (2015) and their linear programming formulation can therefore be reduced from $O(m^2 n^4)$ constraints to $O(m^2 n^3)$.

Table 1 Comparison for joint-constrained cases by the running time in seconds of NLS and CPLEX

Consider the following $n + 1$ lines: $f_0(u) = 0, f_j(u) = v_{ij} (r_{ij} - u)$. Without loss of generality, we assume $v_{i1} \le v_{i2} \le \dots \le v_{in}$; if some $v_{ij} = v_{i(j + 1)}$, then order $v_{ij} r_{ij} \le v_{i(j + 1)} r_{i(j + 1)}$.

There are $q \le \frac{1}{2} n(n + 1)$ crosspoints for those $(n + 1)$ lines, and the corresponding horizontal coordinate is denoted as $I_1 \le I_2 \le \dots \le I_q$. For any consecutive $I_k < I_{k + 1}$ with different horizontal coordinates, the relative order of all lines does not change between $u \in (I_k, I_{k + 1})$ because for any $j', j''$, the relative order $f_{j'}(u)$ and $f_{j''}(u)$ does not change between $u \in (I_k, I_{k + 1})$. Scan from $u \rightarrow -\infty$ to $+\infty$. Initially, when $u \rightarrow -\infty$, the order from highest to lowest must be $(n, n - 1, \dots , 1, 0)$. The following facts provide us the intuition.

• If there is a crosspoint u only for two lines $j_1 > j_2$, i.e., $f_{j_1}(u) = f_{j_2}(u)$, their ranks must be consecutive and should be swapped before and after u, i.e., $f_{j_1}(u - \varepsilon ) > f_{j_2}(u - \varepsilon )$ and $f_{j_1}(u + \varepsilon ) < f_{j_2}(u + \varepsilon )$.

• When $r \ge 2$ lines $j_1> \dots > j_r$ share the same crosspoint u, i.e., $f_{j_1}(u) = \dots = f_{j_r}(u)$, their ranks must be consecutive at $u - \varepsilon$ and should be reversed at $u + \varepsilon$.

Assume $j_1> \dots > j_r$, which implies $f_{j_1}(u - \varepsilon )> \dots > f_{j_r}(u - \varepsilon )$. The common crosspoint can be regarded as $\frac{1}{2}r(r - 1)$ overlapped crosspoints. The reverse operation can be decomposed to a series of $\frac{1}{2}r(r - 1)$ swap operations like the well-known bubble sort,

$$\begin{aligned} (j_1, j_2), (j_1, j_3), \dots , (j_1, j_r), (j_2, j_3), \dots , (j_2, j_r), \dots , (j_{r - 1}, j_r). \end{aligned}$$

Note that each swap only involves two consecutive-ranked lines. After this procedure, the order is reversed, which exactly matches $f_{j_r}(u + \varepsilon )> \dots > f_{j_1}(u + \varepsilon )$.

To build an algorithm, we collect all q crosspoints $(u, j', j'')$ such that $f_{j'}(u) = f_{j''}(u)$ and $j' < j''$ into a set $\mathcal {I} = \{(u_k, j'_k, j''_k)\}_{k = 1}^q$. Sort all crosspoints to satisfy the following condition for all $k = 1, \dots , (q - 1)$:

Rule 1

• $u_k \le u_{k + 1}$;

• If $u_k = u_{k + 1}$, then $j'_k \ge j'_{k + 1}$;

• If $u_k = u_{k + 1}$ and $j'_k = j'_{k + 1}$, then $j''_k > j''_{k + 1}$.

Starting from the order permutation $\{n, \dots , 0\}$ as $u \rightarrow -\infty$, the above discussion ensures that we can follow $k = 1, \dots , q$, and for each step we swap two elements $j'_k$ and $j''_k$ and record the current order permutation. Note the merit is that $j'_k$ and $j''_k$ must be consecutive. After we finish this procedure, we have $q + 1 = O(n^2)$ permutations.

Obviously the crosspoint with smaller u should be processed earlier since we are simulating u from $-\infty$ to $+\infty$. To illustrate how this condition works under crosspoints with same u, Fig. 1 consists of 5 lines $\{1, 2, 3, 4, 5\}$ ranked from lowest absolute slope to highest. There are 4 crosspoints at $u = 0.5$ (3 ones are overlapped). For $u = 0.5 - \varepsilon$, the order is (5, 3, 2, 4, 1); after a series of swaption: (3, 5), (2, 5), (2, 3), (1, 4), the order becomes (2, 3, 5, 1, 4), which is exactly the order at $u = 0.5 + \varepsilon$.

Fig. 1

An example for the ordering condition

Since there are $O(n^2)$ crosspoints, $\vert \mathcal {T}_i(c) \vert = O(n^2)$. However, we can claim a tighter bound for $\sum _{c = 0}^n \vert \mathcal {T}_i(c) \vert$.

Lemma 1

$\sum _{c = 0}^n \vert \mathcal {T}_i(c) \vert = O(n^2)$.

Proof

In our procedure, we swap two consecutive elements in a permutation $O(n^2)$ times. When swapping $j'_k$ and $j''_k$, which are ranked $(c + 1)$ and c respectively, the set of top-c updates while all sets of top-$c'$ remain same for $c' \ne c$. There are three cases:

• $j'_{k} = 0$. We can extract the top-$(c - 1)$ elements to $\mathcal {T}_i(c')$ for all $c' \ge c$. This case can occur exactly n times so at most $O(n^2)$ elements are created.

• The order of $j'_k$ is higher than 0. In this case, a swap will create one element belonging to $\mathcal {T}_i(c)$.

• The order of $j'_k$ is lower than 0. No new elements are created.

Therefore, $O(n^2)$ swaps will finally make $\sum _{c = 0}^n \vert \mathcal {T}_i(c) \vert = O(n^2)$. $\square$

Algorithm 2 is a summary of the above procedure and Theorem 2 analyzes the complexity.

Theorem 2

In Algorithm 2, if we record $V_i(\hat{S}_c)$ and $R_i(\hat{S}_c)$ instead of each assortment $\hat{S}_c$ for each nest i, the algorithm terminates in $O(n^2 \log n)$ time with all candidate sets $\mathcal {T}_i(c)$ for $c = 1, \dots , n$.

Proof

The dominant part of preprocessing is generating and sorting $O(n^2)$ crosspoints which costs $O(n^2 \log n)$. For each crosspoint involving top-c, we swap the corresponding elements in the permutation as well as their ranks. If we use $(v_{ij}, v_{ij} r_{ij})$ instead of each $j = 0, \dots , n, V_i(\hat{S}_c)$ and $R_i(\hat{S}_c)$ are prefix sum queries which can be solved in $O(\log n)$ per query, by some data structures as simple as Fenwick Trees (Fenwick 1994). Therefore, it is sufficient to finish the process in $O(n^2 \log n)$ time. $\square$

Generating all $\mathcal {T}_i(c)$ for $i = 1, \dots , m$ and $c = 1, \dots , n$ can be finished in $O(mn^2 \log n)$.

Fig. 2

The maximization of $g_i^{k - 1}(z)$ and $a^i_k - b^i_k z$

1.2 B Forming the piecewise-linear functions

In this note, an important technical issue is to build the analytical form of piecewise-linear functions $g_i(z) = \max _j \{a^i_j - b^i_j z\}$ and sum them up into $G(z) = \sum _i g_i(z)$.

Formally, $g_i(z)$ is a maximization of a sequence of linear functions $\mathcal {L} = \{a^i_j - b^i_j z\}$, simplified as $g_i(z) := \max _{j = 1}^{\vert \mathcal {L} \vert } \{a^i_j - b^i_j z\}$. Without loss of generality, we assume $b^i_1> \dots > b^i_{\vert \mathcal {L} \vert }$ (ranking the lines from highest steepness to lowest requires an $O(\vert \mathcal {L} \vert \log \vert \mathcal {L} \vert )$ sorting procedure) and no $(a^i_{j'}, b^i_{j'})$ such that $\exists b^i_{j'} = b^i_{j''}$ and $a^i_{j'} < a^i_{j''}$ since line $j'$ is redundant.

With a slight abuse of notation, we define $g_i^k(z) = \max _{j = 1}^k \{a^i_j - b^i_j z\}$. For $k \ge 2$, the recursive formulation is $g_i^k(z) = \max \{g_i^{k - 1}(z), a^i_k - b^i_k z\}$. We know that the (sub)gradient of $g_i^{k - 1}(z)$ for any $z \in {\mathbf {R}}$ is strictly less than $-b^i_k$. There must be a breakpoint u such that $g_i^{k - 1}(u) = a^i_k - b^i_k u, z < u \Leftrightarrow g_i^{k - 1}(z) > a^i_k - b^i_k z$ and $z > u \Leftrightarrow g_i^{k - 1}(z) < a^i_k - b^i_k z$ (see Fig. 2 for an illustration). Therefore $g_i^k(z)$ can be formulated by some leftmost linear piece of $g_i^{k - 1}(\cdot )$ plus a new linear piece

$$\begin{aligned} g_i^k(z) = \left\{ \begin{array}{ll} g_i^{k - 1}(z), &{}\quad z \le u\\ a^i_k - b^i_k z, &{}\quad z \ge u \end{array} \right. \end{aligned}$$

After repeatedly deleting the rightmost linear piece until we find the intersection point u, we create a new linear piece for $z \ge u$ and the function is then updated to $g_i^k(\cdot )$. The iteration ends up with $g_i(z) = g_i^{\vert \mathcal {L} \vert }(z)$.

The final step is to sum $\{g_i(z) := a_i(z) - b_i(z) z\}$ up into the function $G(z) := A(z) - B(z) z$. Obviously, when $z \rightarrow -\infty , G(z) = -v_0 z + \sum _{i = 1}^m (a_i(z) - b_i(z) z)$. So we can start from $A(z) := \sum _{i = 1}^m a_i(z)$ and $B(z) := v_0 + \sum _{i = 1}^m b_i(z)$ where $z \rightarrow -\infty$. The set of $G(\cdot )$’s breakpoints consists of all breakpoints of $g_i(\cdot )$ for all $i = 1, \dots , m$; also $A(\cdot )$ and $B(\cdot )$. With the ascending order of z, each breakpoint of $g_i$ such that $g_i(z - \varepsilon ) := a_i(z - \varepsilon ) - b_i(z - \varepsilon ) z$ and $g_i(z + \varepsilon ) := a_i(z + \varepsilon ) - b_i(z + \varepsilon ) z$ makes a delta change that $A(z) := A(z - \varepsilon ) + a_i(z + \varepsilon ) - a_i(z - \varepsilon )$ and similarly B(z). The procedure only involves linear scanning and sorting so $O(mn^2 \log mn)$ is sufficient to build the piecewise-linear function.