Consider the integers $\{1,\dots, N\}$ for some positive integer $N$. Let us suppose that for each $\{1, \dots, N\}$ there is an associated probability $p_1, \dots, p_N$. We also define an integer threshold $1 \leq n < N$.

We sample independently and repeatedly without replacement from $\{1,\dots, N\}$. For each $i \in \{1,\dots, N\} $ we sample the integer $i$ with probability $p_i$. As we are sampling without replacement these probabilities will have to be appropriately scaled after each new integer is sampled.

Is there a fast algorithm to compute the expected number of distinct integers less than or equal to the threshold $n$ in a sample of size $x$?

The expected value we want to compute is $\sum_{i=1}^n q_i$ where $$q_i = \sum_{y = 1}^x \sum_{\substack{\{j_k \neq i\text{ distinct}\}\\k=1,\dots,y-1}}\left(\prod_{k=1}^{y-1}\frac{p_{j_k}}{1 - \left(\sum_{\ell = 1}^{k-1} p_{j_{\ell}}\right)}\right)\frac{p_i}{1 - \left(\sum_{\ell = 1}^{y-1} p_{j_{\ell}}\right)}.$$

However, this is infeasible to compute for anything but the smallest problem instances.

Is there an efficient algorithm for this problem?

  • $\begingroup$ Do you mean probability proportional to $p_i$ rather than equal to $p_i$? Even if initially we have $\sum_i p_i=1$ after the first sampling the sum of the remaining $p_i$ will not be 1 anymore. $\endgroup$ Apr 17, 2021 at 19:12
  • $\begingroup$ @MaxAlekseyev Yes. After every integer is sampled the probability of the remaining integers will of course have to be scaled to make them add to one. $\endgroup$
    – Simd
    Apr 17, 2021 at 19:33

1 Answer 1


I'm not sure why you ask for "distinct integers" when sampling without replacement guarantees distinctness.

Let $q_i=1-p_i$. The ordinary generating function $$F(u,y) = \prod_{i=1}^n (q_i+p_i uy) \prod_{i=n+1}^N (q_i+p_iy)$$ counts subsets with $u$ monitoring the choices in $[1..n]$ and $y$ monitoring all choices.

To extract the expectation of the number of choices in $[1..n]$, differentiate with respect to $u$. If your sample size is $k$ (sorry but I can't use $x$ for an integer parameter when playing with generating functions), extract the coefficient of $y^k$.

So $$\text{answer} = \frac{[y^k] F_u(1,y)}{[y^k] F(1,y)},$$ where $[y^k]$ means extract the coefficient of $y^k$.

Now to computation. Start with the denominator $$F(1,y) = \prod_{i=1}^N (q_i+p_i y).$$ I don't know of any super-clever way to extract the coefficient of $y^k$ but you can just start with the polynomial $1$ and multiply by one factor at a time discarding any powers greater than $k$. That will have complexity $O(kN)$. Probably there are options with lower theoretical complexity but the constant here is so small that it will be hard to beat for practical sizes.

The numerator is $$F_u(1,y) = \sum_{i=1}^n p_i y\prod_{j\ne i} (q_j+p_jy),$$ which you can compute in similar manner.

ADDED. I'll explain how to get a good estimate for large problems, assuming that $k$ is not very small. I'll consider only the denominator. The function $F(1,y)$ is the pgf of the sum of $N$ independent Bernoulli random variables. The sum satisfies a local central limit theorem (LCLT) but that won't help immediately unless $k$ is close to the mean of the sum. Now introduce a parameter $\alpha\gt 0$ and replace $q_i+p_iy$ by $(q_i + \alpha p_iy)/(q_i+\alpha p_i)$ for each $i$. Solve for $\alpha$: $$ \sum_{i=1}^N \frac{\alpha p_i}{q_i+\alpha p_i}=k,$$ which should be fairly easy numerically as there is no cancellation and the left side is monotonically increasing in $\alpha$.

With this value of $\alpha$, the mean of the sum of the tilted Bernoulli random variables with pgf $(q_i + \alpha p_iy)/(q_i+\alpha p_i)$ is exactly $k$. Now you can apply a LCLT such as Theorem 6.3 in this paper to get quite a good estimate.

If you want an arbitrarily accurate value, you can apply a fourier inversion theorem to the characteristic function and evaluate it numerically to the precision you desire. I don't know if this will be easy, but for sure it will be advantageous to tilt the distribution first as I did above.

  • $\begingroup$ I'm not sure why you ask for "distinct integers" ... That surely comes from a previous (and equivalent) problem statement, where one samples with replacements (leaving the probabilities unchanged) and stops after getting $n$ different values. math.stackexchange.com/questions/4096509/… $\endgroup$
    – leonbloy
    May 11, 2021 at 1:45
  • $\begingroup$ Yes probably, but "without replacement" and "distinct" are the same, which could have indicated a mistake in the description. $\endgroup$ May 11, 2021 at 6:09
  • $\begingroup$ @leonbloy Yes it was exactly that. $\endgroup$
    – Simd
    May 11, 2021 at 15:09
  • $\begingroup$ @BrendanMcKay thank you for this answer. I would be very interested in any methods with lower theoretical complexity. In particular that might be feasible when $N \approx 10^8$ and $k \approx 10^5$. $\endgroup$
    – Simd
    May 11, 2021 at 15:11
  • $\begingroup$ Thank you again. $\endgroup$
    – Simd
    May 22, 2021 at 12:05

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.