Image detection using combinatorial auction

Consider an auction of a variety of different goods. Unlike a classical auction, where goods are sold individually, a bidder can place a bid on a bundle of goods in a combinatorial auction. The auction manager receives multiple price offers for different bundles of goods. The goal of the auction manager is to maximize his revenue by allocating the goods to the highest submitted bids. Naturally, the same good may not be allocated to multiple bidders. This means that if two bidders place bids over bundles with overlapping goods, only one of the bids can win.

Let $\mathcal{G}=\{\gamma_{1},\ldots,\gamma_{A}\}$ be a set of given goods, and let $\mathcal{B}=\{\beta_{1},\ldots,\beta_{B}\}$ be a set of bids. Each bid $\beta_{i}$ includes a bundle of goods $g(\beta_{i})\subseteq\mathcal{G}$, which the bidder seeks to acquire. For this bundle, the bidder of $\beta_{i}$ proposes a price, $p(\beta_{i})$. Let $x_{i}$ be an indicator variable that corresponds to the winning bids. That is, we set $x_{i}=1$ if $\beta_{i}$ is a winning bid and $x_{i}=0$ otherwise. An allocation, $\bm{\pi}$, is a set of winning bids that do not overlap. Specifically, $\bm{\pi}\subseteq\mathcal{B}$ is a subset of the bids, where $\forall\beta_{2}\neq\beta_{1}\in\bm{\pi},g(\beta_{1})\cap g(\beta_{2})=\varnothing$. Therefore, an allocation cannot include two bids that include the same good.

Using this notation, the combinatorial auction problem is given by:

The first constraint ensures that each good can be allocated once at most, while the second constraint guarantees that each bid can be either included or excluded from the chosen allocation. We define the optimal allocation as the allocation that maximizes the revenue of the auction manager. The problem of determining the optimal allocation is known as the WDP [23]. It is known that the WDP is an NP-hard problem and cannot be approximated, within any constant, by a polynomial-time algorithm [32, 33].

We now cast the constrained optimization problem (4) as a variation of the WDP. We associate every pixel in the measurement with a single good and let $\mathcal{B}$ be a set of bids. Accepting $\beta_{ij}\in\mathcal{B}$ implies that our estimator includes an image occurrence whose upper left corner is located at the $(i,j)$-th pixel. The goods requested by the bid, $g(\beta_{ij})$, are a $W\times W$ set of pixels, corresponding to a possible image occurrence in the measurement. The offered price of the bid, $p(\beta_{ij})$, is the correlation between the measurement and the template image, namely,

The revenue of the auctioneer is the sum of offered prices by all accepted bids. Since we are estimating the upper left corner coordinates of images in the measurement, we require that the estimated coordinates will satisfy $\hat{c}^{(k)}_{n}\leq N-W+1$ and $\hat{c}^{(k)}_{m}\leq M-W+1$, so we do not exceed the measurement’s dimension. We denote $\tilde{N}=N-W+1$ and $\tilde{M}=M-W+1$.

Let us define $\mathbf{X}\in\{0,1\}^{\tilde{N}\times\tilde{M}}$ as a matrix of indicators, in which $x_{ij}\in\mathbf{X}$ is equal to $1$ if the bid $\beta_{ij}$ is included in the allocation, and $0$ otherwise. To enforce the separation constraint (2), we require that

where $\mathbbm{1}_{W}$ is an all-ones $W\times 1$ column vector and $\tilde{\mathbf{X}}_{ij}$ is a $W\times W$ sub-matrix of $\mathbf{X}$, whose upper left entry is the $(i,j)$ entry in the matrix $\mathbf{X}$. Using this notation, and assuming $K$ is known, the constrained likelihood problem (4) can be written as:

where the last two conditions are, together, equivalent to (7). Therefore, solving (8) is equivalent to (4), where the non-zero entries of the optimal $\mathbf{X}$ correspond to the MLE of the locations in (4).

We underscore that the optimization problem (8) is slightly different from the classical combinatorial auction (5). In the classical setup, the auctioneer aims to include as many bids as possible in the optimal allocation to maximize the revenue, whereas in our model the optimal allocation must include exactly $K$ bids. Furthermore, in our model, $p(\beta_{ij})$ can be negative. Thus, we cannot use existing algorithms that solve the classical WDP (5), such as the one proposed in [23], and we design a new algorithm tailored for our constrained setup. Importantly, the designed algorithm must be computationally efficient: solving (8) in a brute-force manner requires searching over an exponential number of bids, rendering this approach intractable. For example, if $N,M=30$ and $K=5$, then the number of possible allocations is $\sim 10^{12}$. The next section introduces a pruning technique to solve (8), which dramatically reduces the number of explored allocations.

The algorithm begins with forming a binary tree data structure. Binary trees, widely used in search algorithms in computer science [2], are suitable for examining every possible allocation through a systematic navigation to determine the optimal allocation.

Given a set of bids, each node in the tree corresponds to a simple binary question: whether to include or exclude the bid in the allocation. This results in a full binary tree where, except for the leaf nodes, each node has two children corresponding to the same bid. We introduce a simple example to illustrate the construction of the binary tree. Let us consider $y\in\mathbb{R}^{4}$, a template image of a unit size, $s\in\mathbb{R}$, and $K=2$. As explained in Section 2, we cast the image detection problem as a variation of the WDP. In this example, every bid includes a bundle of a single good (since the template image is of a unit size). We first set $\beta_{1}$ as the tree’s root, and $\beta_{4}$ forms the tree’s leaves. To explore an allocation that includes the current node, we navigate to the node’s left child and add the current node to our partial allocation, $\bm{\pi}_{k}$. Conversely, navigating to the right child represents the alternative choice of excluding the current node from $\bm{\pi}_{k}$. Figure 3 demonstrates the binary tree we obtain in this simple example.

As previously mentioned, the binary tree enables us to conduct an exhaustive search across all possible allocations. Therefore, by using this search mechanism, the optimal allocation will inevitably be encountered during our search.

Exhaustively searching through the tree leads necessarily to the optimal allocation, but is computationally infeasible even for problems of small dimensions. To overcome this issue, we use the algorithmic technique of branch-and-bound and present a pruning method to avoid exploring a large portion of the allocations in advance by eliminating sub-trees that cannot yield the optimal allocation.

We illustrate this method using the example of Figure 3. We refer to a smaller binary tree, spanned by a node other than the root, as a sub-tree. Let us assume that our optimal allocation is $\bm{\pi}_{\text{opt}}=\{\beta_{1},\beta_{3}\}$ and that we have already encountered our optimal allocation $\bm{\pi}_{\text{opt}}$ in our search. At this point, we are unaware that this is the optimal allocation, so we continue exploring the remaining possible allocations. Assume we have chosen to exclude $\beta_{1}$ from the empty partial allocation, $\bm{\pi}_{0}$. Following that, two bids are left to be found in the sub-tree spanned by $\beta_{2}$. Figure 4 demonstrates this sub-tree. The maximal revenue achievable by exploring this sub-tree is bounded by the sum of the prices of the two bids with the highest revenue within it. Since this revenue is less than what is achieved by the optimal allocation (by the assumption that we already encountered the optimal allocation), we can avoid exploring the allocations in this sub-tree. Thus, in this example, we reduce the computational burden by half, as we prune half of the binary tree.

We now formulate this example as a general pruning mechanism. Say we have reached a partial allocation $\bm{\pi}_{k}$. We sort the bids in $\overline{C}_{\bm{\pi}_{k}}$ (the set of allowed bids) in descending order according to their revenues. We denote the $K-k$ maximal revenue bids in $\overline{C}_{\bm{\pi}_{k}}$ as $\overline{C}^{(K-k)}_{\bm{\pi}_{k}}$. Summing the revenues of the bids in $\overline{C}^{(K-k)}_{\bm{\pi}_{k}}$ bounds on the remaining revenue to be allocated in this sub-tree

This function will enable us to determine if the optimal allocation cannot be achieved within the current sub-tree, and thus we can prune it. Specifically, we prune the sub-tree if the sum of $p(\bm{\pi}_{k})$ and $h(\bm{\pi}_{k})$ is less than the revenue of our current optimal allocation, $p(\bm{\pi}_{\text{opt}})$, i.e.,

Since we only prune sub-trees that do not lead to the optimal allocation, our proposed algorithm remains optimal. We note, however, that $h(\bm{\pi}_{k})$ is not necessarily a tight upper bound as it ignores the constraint that the remaining $K-k$ bids must not conflict with each other.

Our algorithm examines each potential optimal allocation, and thus the optimality of the algorithm is unaffected by rearranging the order of the bids. We aim to reorder the bids in such a way that maximizes the efficiency of the pruning mechanism. Note that the earlier we secure a high revenue allocation $p(\tilde{\bm{\pi}}_{\text{opt}})$, the more the stopping criterion (10) is met, resulting in fewer allocations to be examined. Hence, we begin the algorithm by sorting the bids by their revenue in descending order. The bid with the highest revenue is fixed as the tree’s root, while the bid with the lowest revenue forms the tree’s leaves. This strategy enhances the frequency in which the pruning condition is met, thereby boosting its efficiency and accelerating the algorithm’s runtime. For example, in the experiments in Section 4, we observed that sorting the bids can lead to acceleration by a factor of around 10.

We are now ready to describe our algorithm for detecting image occurrences in a noisy measurement. We will use the example above to explain the guiding principles. We start by sorting the bids in descending order by their revenues, as described in Section 3.3. An example of a sorted tree is presented in Figure 5.

Our algorithm is designed such that when we encounter a node, it is added to our partial allocation only if it does not conflict with the partial allocation and if the partial allocation itself does not satisfy the pruning condition (10). We begin our search at the tree’s root, $\beta_{1}$. As this is the first node we encounter, we add the node to our partial allocation and move to the left child, $\beta_{3}$. At this point, $\beta_{3}$ is added to the partial allocation, since it does not conflict with our partial allocation nor does the partial allocation itself satisfy the pruning condition (10). Since this is the first allocation we attain, this is the current optimal allocation. We note that the first allocation we explore is the output of the greedy algorithm. This means that the greedy algorithm is not necessarily optimal, since the optimal allocation may be achieved in a later stage. In our specific simple example, indeed the greedy algorithm attains the optimal solution (since the first two bids do not conflict), but in order to explain the full mechanism of the algorithm, we describe the remaining steps.

We then remove the last bid we added, $\beta_{3}$, and continue examining allocations excluding that bid. To do so, we move to the right child of $\beta_{3}$. Upon reaching $\beta_{2}$, we check if our partial allocation $\bm{\pi}_{1}$ satisfies the pruning condition (10). Since the optimal allocation is $\bm{\pi}_{\text{opt}}=\{\beta_{1},\beta_{3}\}$, the pruning condition is met, and hence we prune the sub-tree spanned by $\beta_{2}$ and return to the parent, $\beta_{3}$. Since we examined the sub-tree spanned by $\beta_{3}$, we next return to the parent, $\beta_{1}$, and remove $\beta_{1}$ from our partial allocation. This way, we examined all allocations that include $\beta_{1}$, and we move to the node’s right child. This process continues until the algorithm explores the entire tree. Upon completion, it returns the optimal allocation. The general algorithmic pipeline is described in Algorithm 1.

In many real-world scenarios, such as electron microscopy [1, 11], the exact number of image occurrences $K$ is unknown. This imposes a major challenge given that both the greedy algorithm and Algorithm 1 require a fixed $K$. A common practice is solving (4) for various $K$ values, choosing the one that induces the steepest change in the objective value. This heuristic is called the “knee” method and has been widely adopted across diverse domains, such as clustering [37] and regularization [18].

Here, we propose employing the gap statistics principle, which can be thought of as a statistical technique to identify the “knee” [37]. This principle is widely applied to a variety of domains, as discussed for example in [26, 27]. The underlying idea of gap statistics is to adjust the objective value curve (as a function of the possible values of $K$) by comparing it with its expectation under a null reference. Let us define the gap statistic as:

where $p(\bm{\pi}_{\text{opt}}(K))$ represents the revenue of the optimal allocation for $K$ image occurrences (computed by Algorithm 1), and $\mathbb{E}^{*}p(\bm{\pi}_{\text{opt}}(K))$ is the expectation for a measurement of dimensions $N\times M$ derived from a “null” reference distribution. To form the “null” distribution, we randomly rearrange the pixels of the measurement, resulting in an unstructured measurement, with no template image occurrences. Specifically, we estimate the expectation under the null by $E^{*}p(\bm{\pi}_{\text{opt}}(K))\approx\frac{1}{R}\sum_{r=1}^{R}p(\bm{\pi}_{\text{opt}}^{r}(K))$, where $p(\bm{\pi}_{\text{opt}}^{r}(K))$ is the revenue achieved by the optimal allocation, at the $r$-th rearrangement; in the experiments below, we set $R=50$. The estimate of $K$—the number of image occurrences—is the $K$ that maximizes $\text{gap}(K)$ and is denoted by $\hat{K}$.

In the first experiment, we assume that the number of image occurrences $K$ is known, and both algorithms receive the true number of image occurrences $K$ as input. We set $N=M=40$ and $K=4$, and choose a template image with all-ones entries of size $W=3$. We generate the measurement such that the separation condition (2) is met, and at least one pair of image occurrences is separated precisely by $W$ pixels. For each SNR level, we conduct 1000 trials, each with a fresh measurement. Figure 6 shows the average $F_{1}$ score of Algorithm 1 and the greedy algorithm as a function of the SNR. Evidently, Algorithm 1 outperforms the greedy algorithm in all SNR regimes. As expected, in high SNR regimes, Algorithm 1 achieves $\text{F}_{1}=1$. On the contrary, the greedy algorithm fails even in the high SNR regime, achieving $\text{F}_{1}=0.88$, due to the proximity of the image occurrences in the measurement. The average runtime per trial of Algorithm 1 was approximately 80 seconds, whereas the average runtime of the greedy algorithm was 0.1 seconds.

We repeated the same experiment, but now the image occurrences are well-separated according to (15). The results of this experiment are illustrated in Figure 7. As expected, in this case, both algorithms achieve the same results, and the running time remains the same.

We repeated the first experiment from the previous section (when the image occurrences satisfy (2) but not (15)), but now the number of image occurrences $K$ is unknown. To estimate $K$, we used the gap statistics principle. This estimate is then used as an input for both algorithms to determine the image locations. The results are presented in Figure 8.

Similar to the previous case, in the unknown $K$ scenario, Algorithm 1 outperforms the greedy algorithm in all SNR regimes. In addition, Algorithm 1 provides better estimates of the number of signal occurrences $K$. For example, for SNR$=-12.5$, Algorithm 1 gets $87.2\%$ accuracy in estimating $K$ precisely, while the greedy algorithm gets only $80\%$. This gap drops with the SNR until they eventually converge to the same estimations when the SNR is low enough.

Our first experiment is conducted on the EMPIAR-10081 data set of the Human HCN1 Hyperpolarization-Activated Channel [22]. The measurement is presented in Figure 9. We set the radius of the disc used as a template image to 3 pixels. The number of image occurrences is estimated based on gap statistics. Figure 10 shows the output of both algorithms. Algorithm 1 identifies two separated images. In contrast, the greedy algorithm identifies a single occurrence. For comparison, we also show the output of a popular detection algorithm for cryo-EM data sets (particle picker), called TOPAZ [3]. This algorithm is based on a deep learning technique (see also [7]).

Our second experiment is conducted on the EMPIAR-10217 data set of the bovine liver glutamate dehydrogenase [21]. We set the radius of the disc used as a template image to 6 pixels. The measurement is presented in Figure 9 and the results of both algorithms are displayed in Figure 11. Similar to the previous example, Algorithm 1 identifies two different particle images, while the greedy algorithm finds only one particle image in the measurement.