How to (efficiently) compute the solution to the following equation?

34 Views Asked by At

\begin{align}\sum_{i=1}^N\min(1, w^{(i)}/\alpha ) = M,\end{align} where $\alpha \in \mathbb{R}^{+}$ is strictly positive, unknown and has to be found, the $w^{(1)}, \dots w^{(N)}$ are positive, known and sum to $1$ and $M \in \mathbb{N}$ is known, too.

This is a step in an algorithm called the Stratified Optimal Resampling (SOR) algorithm due to Fearnhead and Liu (2007), see here: eprints.lancs.ac.uk/745/. Unfortunately, they do not mention how to solve this equation for $\alpha$, but I want to solve this problem in as little computing time as possible. I don't really see an easy way to do this, though. What am I overlooking here?

1

There are 1 best solutions below

3
On BEST ANSWER

The left hand side is piecewise linear in $\alpha^{-1}.$ You can efficiently find a solution on each piece if one exists.

The pieces can be traversed efficiently given access to the sorted $w^{(i)}.$ I’m assuming your data isn’t so big that sorting is an issue.

You can make this more efficient by keeping track of partial sums; note the left hand side can be expressed in terms of $\alpha$ and the sum of elements less than $\alpha.$

The left hand side is also increasing in $\alpha^{-1}$ which gives another way to find a solution: binary chop.

Using piecewise linearity

Assume $w^{(1)}\leq \dots \leq w^{(N)}.$ For $1\leq k\leq n-1$ define $$ f_k(\alpha) = N-k+\frac1\alpha\sum_{i=1}^k w^{(i)}.$$

Then $f(\alpha)=f_k(\alpha)$ whenever $w^{(k)}\leq \alpha\leq w^{(k+1)}.$ The unique solution $\alpha_k$ of $f_k(\alpha_k)=M$ is $$\alpha_k=\left(\sum_{i=1}^k w^{(i)}\right)\Bigg/(M-N+k).$$ If $\alpha_k\in[w^{(k)},w^{(k+1)}],$ then $f(\alpha_k)=M.$ Conversely if $f(w^{(1)})<M<f(w^{(N)})$ then the solution of $f(\alpha)=M$ lies in some interval $[w^{(k)},w^{(k+1)}],$ so can be found by checking each $\alpha_k.$ This might be implemented in practice as a loop for $k$ from $1$ to $n-1,$ accumulating $\sum_{i=1}^k w^{(i)}.$

Using binary chop

A less efficient but perhaps more readable solution is to implement $f(x)$ in the simplest possible way (should be linear time in $N$), then perform a binary search between some bounds ($\min_i w^{(i)}$ and $\max_i w^{(i)}$?) for the value of $\alpha$ that makes $f(\alpha)=M.$ This will get you $d$ binary digits in time $O(Nd).$