# Combinatorics - calculating probability of choosing n people from m groups of k people

I am trying to calculate and understand the following question: I have 10 groups of 10 people each and am choosing 10 people at random. What is the probability that all 10 people I have chosen are from one group? From 2 groups? Etc.

• Let's try a simpler example: $3$ groups of $3$ people and you choose $3$ of the $3^2=9$ people. There are ${9 \choose 3}=84$ equally likely ways of choosing the $3$ in total. Can you see why there are $3$ ways of choosing them from a single group, $54$ from two groups and $27$ ways of choosing them from three groups? Commented Feb 8, 2022 at 15:41

Consider $$m=10$$ groups of size $$n=10$$ each. Obtain a sample of size $$N$$ by withdrawing one element at a time. For the possible counts $$0, 1, 2, \ldots, m,$$ let $$p_N(i)$$ be the chance that exactly $$i$$ groups are represented in this sample.

Notice, now, that when $$i$$ groups are represented in a sample of size $$N,$$ there remain $$nm-N$$ unsampled elements, of which $$in-N$$ are in the groups appearing in the sample. The chance of selecting one of these at the next step is therefore $$(in-N)/(nm-N)$$, while the chance of selecting a different one--and thereby increasing $$i$$ to $$i+1$$--is $$1 - (in-N)/(nm-N).$$ Consequently, the $$p_N(i)$$ enjoy the recursive relationship

$$p_{N+1}(i) = p_N(i) \left(\frac{ni-N}{nm-N}\right) + p_N(i-1)\left(1 - \frac{n(i-1)-N}{nm-N}\right).$$

At the beginning, with a sample size of $$1,$$ clearly $$p_1(1)=1$$ and all the other $$p_1(i)=0.$$ Updating this vector of length $$m+1$$ $$N-1$$ times requires only $$(m+1)(N-1)$$ of the foregoing calculations. That's efficient. (Well before $$m$$ and $$n$$ become large enough to make this onerous, asymptotic approximations will become highly accurate.)

Here is R code to carry out this simple calculation:

n <- 10 # Group size
m <- 10 # Number of groups
N <- 10 # Sample size

P <- c(0, 1, rep.int(0, m-1))
i <- seq_len(m)
for (k in seq_len(N)[-1] - 1)
P <- c(0, (P[i + 1] * (n * i - k) + P[i] * n * (m - i + 1)) / (n * m - k))


As a test, let's compare it to a direct, brute-force simulation of a half million samples.

set.seed(17)
pop <- rep.int(seq_len(m), n)                             # The elements to sample
x <- replicate(5e5, sum(tabulate(sample(pop, N), m) > 0)) # Take the samples
simulation <- c(0, tabulate(x, m) / length(x))            # Compute frequencies
rbind(simulation, P)                                      # Compare to the calculation

           0        1       2        3       4      5     6     7     8      9       10
simulation 0 0.00e+00 0.0e+00 0.000226 0.00865 0.0928 0.313 0.389 0.172 0.0234 0.000572
P          0 5.78e-13 4.8e-07 0.000204 0.00884 0.0922 0.314 0.389 0.172 0.0234 0.000578


The headers are the group counts. You can see the agreement is good.

This is really a hypergeometric problem. If you have 10 groups with 10 people in each, and you draw randomly without replacement 10 people, let $$X_i$$ be the number of people drawn from group $$i$$. Then the probability distribution is given by $$\DeclareMathOperator{\P}{\mathbb{P}} \P\left\{ X_i=x_i, i=1,2,\dotsc,10 \right\} = \frac{\prod_i^{10} \binom{10}{x_i}}{\binom{100}{10}}$$ but your interest is in the random variable "number of groups" $$G=\sum_i \{X_i > 0 \}$$ where the curly braces is the< indicator function. For $$g$$ small or close to 10, this is a small sum of hypergeometric probabilities, while for <intermediate values it becomes complicated ...

The comment by @Henry gives a numerical exmple with smaller numbers. Let us do it with your group sizes, with the help< of R. For the intermediate values of $$g$$ it is most practical to resort to simulation.

library(extraDistr)
### A helper function to save typing:
myhyper <- function(x) {
p <- length(x)
n <- sum(x)
}

### Probability of exactly one group, G=1:
10*myhyper(c(10, rep(0, 9)))
[1] 5.776904e-13
# We multiply by 10 because there are 10 ways of
# choosing that group
### Probability of G=2

prob2 <- 0.0
for (i in seq(from=1,  to=9,  by=1)) prob2 <- prob2  +
myhyper(c(i, 10-i, rep(0, 8)))

( prob2 <- choose(10, 2) * prob2 )
4.802878e-07

### Then we resort to simulation:

set.seed(7*11*13) # My public seed

prob2_by_sim <- rmvhyper(1E7, rep(10, 10), 10)

### There is one sample in each row. Let us calculate ### the number of zeros, by row:

zeros <- apply(prob2_by_sim, 1,
FUN = function(row) sum(row==0) )

tzeros <- table(zeros)

tzeros/1E7
zeros
0         1         2         3         4         5         6         7
0.0005821 0.0233285 0.1723967 0.3886847 0.3139600 0.0920672 0.0087765 0.0002039
8
0.0000004


Note the excellent correspondence of the simulation for $$G=2$$ (that is, 8 zeros) with the theoretical calculation above.

Note the similarity of this problem with the coupon collector problem, with number of coupons corresponding to number of groups. But this is a hypergeometric version of the coupon collector problem. There is a paper dedicated to this: The Hypergeometric Coupon Collection Problem and its Dual by Sheldon M. Ross.

• I find it confusing that the software's output is labeled by $10-G$ rather than $G.$
– whuber
Commented Apr 25, 2022 at 22:38