PrevNext
Not Frequent
 0/20

Combinatorics

Authors: Jesse Choe, Aadit Ambadkar, Dustin Miao, Mihnea Brebenel, Peng Bai

How to count.

If you've never encountered any combinatorics before, AoPS is a good place to start.

Resources
AoPS

practice problems; set focus to Counting and Probability for module topics

AoPS

good book

Resources

Resources
CPH

module is based off of this

cp-algo
HE

teaches fundamental combinatorics with a practice problem at the end

AoPS

teaches basic combinatorics concepts

AoPS

teaches more advanced combinatorics concepts

CF

a good blog about the inclusion-exclusion principle

If you prefer watching videos instead, here are some options:

Resources
YouTube

playlist by mathemaniac

YouTube

lectures 16-23

YouTube

Errichto video regarding expected value and sums of subsets

Binomial Coefficients

Focus Problem – try your best to solve this problem before continuing!

The binomial coefficient (nk)\binom{n}{k} (pronounced as "nn choose kk" or sometimes written as nCk{}_nC_k) represents the number of ways to choose a subset of kk elements from a set of nn elements. For example, (42)=6\binom{4}{2} = 6, because the set {1,2,3,4}\{1,2,3,4\} has 66 subsets of 22 elements:

{1,2},{1,3},{1,4},{2,3},{2,4},{3,4}\{1, 2\}, \{1, 3\}, \{1, 4\}, \{2, 3\}, \{2, 4\}, \{3, 4\}

There are two ways to calculate binomial coefficients:

Method 1: Pascal's Triangle (Dynamic Programming) - O(n2)\mathcal{O}(n^2)

Binomial coefficients can be recursively calculated as follows:

(nk)=(n1k1)+(n1k)\binom{n}{k} = \binom{n - 1}{k - 1} + \binom{n - 1}{k}

The intuition behind this is to fix an element xx in the set and choose k1k − 1 elements from n1n − 1 elements if xx is included in the set or choose kk elements from n1n − 1 elements, otherwise.

The base cases for the recursion are:

(n0)=(nn)=1\binom{n}{0} = \binom{n}{n} = 1

because there is always exactly one way to construct an empty subset and a subset that contains all the elements.

This recursive formula is commonly known as Pascal's Triangle.

A naive implementation of this would use a recursive formula, like below:

C++

/** @return nCk mod p using naive recursion */
int binomial(int n, int k, int p) {
if (k == 0 || k == n) { return 1; }
return (binomial(n - 1, k - 1, p) + binomial(n - 1, k, p)) % p;
}

Java

/** @return nCk mod p using naive recursion */
public static int binomial(int n, int k, int p) {
if (k == 0 || k == n) { return 1; }
return (binomial(n - 1, k - 1, p) + binomial(n - 1, k, p)) % p;
}

Python

def binomial(n: int, k: int, p: int) -> int:
""":return: nCk mod p using naive recursion"""
if k == 0 or k == n:
return 1
return (binomial(n - 1, k - 1, p) + binomial(n - 1, k, p)) % p

Additionally, we can optimize this from O(2n)\mathcal{O}(2^n) to O(n2)\mathcal{O}(n^2) using dynamic programming (DP) by caching the values of smaller binomials to prevent recalculating the same values over and over again. The code below shows a bottom-up implementation of this.

C++

/** @return nCk mod p using dynamic programming */
int binomial(int n, int k, int p) {
// dp[i][j] stores iCj
vector<vector<int>> dp(n + 1, vector<int>(k + 1, 0));
// base cases described above
for (int i = 0; i <= n; i++) {
/*
* i choose 0 is always 1 since there is exactly one way
* to choose 0 elements from a set of i elements

Java

/** @return nCk mod p using dynamic programming */
public static int binomial(int n, int k, int p) {
// dp[i][j] stores iCj
int[][] dp = new int[n + 1][k + 1];
// base cases described above
for (int i = 0; i <= n; i++) {
/*
* i choose 0 is always 1 since there is exactly one way
* to choose 0 elements from a set of i elements

Python

def binomial(n: int, k, p):
""":return: nCk mod p using dynamic programming"""
# dp[i][j] stores iCj
dp = [[0] * (k + 1) for _ in range(n + 1)]
# base cases described above
for i in range(n + 1):
"""
i choose 0 is always 1 since there is exactly one way
to choose 0 elements from a set of i elements

Method 2: Factorial Definition (Modular Inverses) - O(n+logMOD)\mathcal{O}(n + \log MOD)

Define n!n! as n×(n1)×(n2)×1n \times (n - 1) \times (n - 2) \times \ldots 1. n!n! represents the number of permutations of a set of nn elements. See this AoPS Article for more details.

Another way to calculate binomial coefficients is as follows:

(nk)=n!k!(nk)!\binom{n}{k} = \frac{n!}{k!(n-k)!}

Recall that (nk)\binom{n}{k} also represents the number of ways to choose kk elements from a set of nn elements. One strategy to get all such combinations is to go through all possible permutations of the nn elements, and only pick the first kk elements out of each permutation. There are n!n! ways to do so. However, note the the order of the elements inside and outside the subset does not matter, so the result is divided by k!k! and (nk)!(n − k)!.

Since these binomial coefficients are large, problems typically require us to output the answer modulo a large prime pp such as 109+710^9 + 7. Fortunately, we can use modular inverses to divide n!n! by k!k! and (nk)!(n - k)! modulo pp for any prime pp. Computing inverse factorials online can be very time-consuming. Instead, we can precompute all factorials in O(n)\mathcal{O}(n) time and inverse factorials in O(n+logMOD)\mathcal{O}(n + \log MOD). First, we compute the modular inverse of the largest factorial using binary exponentiation. For the rest, we use the fact that (n!)1(n!)1×(n+1)1×(n+1)((n+1)!)1×(n+1)(n!)^{-1} \equiv (n!)^{-1}\times (n+1)^{-1} \times (n+1) \equiv ((n+1)!)^{-1}\times (n+1). See the code below for the implementation.

C++

const int MAXN = 1e6;
long long fac[MAXN + 1];
long long inv[MAXN + 1];
/** @return x^n modulo m in O(log p) time. */
long long exp(long long x, long long n, long long m) {
x %= m; // note: m * m must be less than 2^63 to avoid ll overflow
long long res = 1;
while (n > 0) {

Java

import java.util.*;
public class BinomialCoefficients {
private static final int MAXN = (int)1e6;
private static long[] fac = new long[MAXN + 1];
private static long[] inv = new long[MAXN + 1];
/** @return x^n modulo m in O(log p) time. */
private static long exp(long x, long n, long m) {
x %= m; // note: m * m must be less than 2^63 to avoid ll overflow

Python

MAXN = 10**6
fac = [0] * (MAXN + 1)
inv = [0] * (MAXN + 1)
def exp(x: int, n: int, m: int) -> int:
""":return: x^n modulo m in O(log p) time."""
x %= m # note: m * m must be less than 2^63 to avoid ll overflow
res = 1

Solution - Binomial Coefficients

The first method for calculating binomial factorials is too slow for this problem since the constraints on aa and bb are (1ba106)(1 \leq b \leq a \leq 10^6) (recall that the first implementation runs in O(n2)\mathcal{O}(n^2) time complexity). However, we can use the second method to answer each of the nn queries in constant time by precomputing factorials and their modular inverses.

C++

#include <iostream>
using namespace std;
using ll = long long;
const int MAXN = 1e6;
const int MOD = 1e9 + 7;
ll fac[MAXN + 1];
ll inv[MAXN + 1];

Java

import java.io.*;
import java.util.*;
public class BinomialCoeffs {
private static final int MAXN = (int)1e6;
private static final int MOD = (int)1e9 + 7;
private static long[] fac = new long[MAXN + 1];
private static long[] inv = new long[MAXN + 1];
public static void main(String[] args) {

Python

MAXN = 10**6
MOD = 10**9 + 7
fac = [0] * (MAXN + 1)
inv = [0] * (MAXN + 1)
Code Snippet: Counting Functions (Click to expand)

Derangements

Focus Problem – try your best to solve this problem before continuing!

The number of derangements of nn numbers, expressed as !n!n, is the number of permutations such that no element appears in its original position. Informally, it is the number of ways nn hats can be returned to nn people such that no person recieves their own hat.

Method 1: Principle of Inclusion-Exclusion

Suppose we had events E1,E2,,EnE_1, E_2, \dots, E_n, where event EiE_i corresponds to person ii recieving their own hat. We would like to calculate n!E1E2Enn! - \lvert E_1 \cup E_2 \cup \dots \cup E_n \rvert.

We subtract from n!n! the number of ways for each event to occur; that is, consider the quantity n!E1E2Enn! - \lvert E_1 \rvert - \lvert E_2 \rvert - \dots - \lvert E_n \rvert. This undercounts, as we are subtracting cases where more than one event occurs too many times. Specifically, for a permutation where at least two events occur, we undercount by one. Thus, add back the number of ways for two events to occur. We can continue this process for every size of subsets of indices. The expression is now of the form:

n!E1E2En=k=1n(1)k(number of permutations with k fixed points)n! - \lvert E_1 \cup E_2 \cup \dots \cup E_n \rvert = \sum_{k = 1}^n (-1)^k \cdot (\text{number of permutations with $k$ fixed points})

For a set size of kk, the number of permutations with at least kk indicies can be computed by choosing a set of size kk that are fixed, and permuting the other indices. In mathematical terms:

(nk)(nk)!=n!k!(nk)!(nk)!=n!k!{n \choose k}(n-k)! = \frac{n!}{k!(n-k)!}(n-k)! = \frac{n!}{k!}

Thus, the problem now becomes computing

n!k=0n(1)kk!n!\sum_{k=0}^n\frac{(-1)^k}{k!}

which can be done in linear time.

C++

#include <atcoder/modint>
#include <bits/stdc++.h>
using namespace std;
using mint = atcoder::modint;
int main() {
int n, m;
cin >> n >> m;
mint::set_mod(m);

Java

import java.util.Scanner;
public class Main {
public static void main(String[] args) {
Scanner scanner = new Scanner(System.in);
int n = scanner.nextInt();
int m = scanner.nextInt();
long c = 1;
for (int i = 1; i <= n; i++) {

Python

n, m = map(int, input().split())
c = 1
for i in range(1, n + 1):
c = (c * i) + (-1 if i % 2 == 1 else 1)
c %= m
c += m
c %= m
print(c, end=" ")
print()

Method 2: Dynamic Programming

Suppose person 1 recieved person ii's hat. There are two cases:

  1. If person ii recieves person 1's hat, then the problem is reduced to a subproblem of size n2n - 2. There are n1n - 1 possibilities for ii in this case, so we add to the current answer (n1)!(n2)(n - 1)\cdot !(n - 2).
  2. If person ii does not recieve person 1's hat, then we can reassign person 1's hat to be person ii's hat (if they recieved person 1's hat, then this would become first case). Thus, this becomes a subproblems with size n1n - 1, are there n1n - 1 ways to choose ii.

Thus, we have

!n=(n1)(!(n2)+!(n1))!n = (n - 1)(!(n - 2) + !(n - 1))

which can be computed in linear time with DP. The base cases are that !0=1!0 = 1 and !1=0!1 = 0.

C++

#include <atcoder/modint>
#include <bits/stdc++.h>
using namespace std;
using mint = atcoder::modint;
int main() {
int n, m;
cin >> n >> m;
mint::set_mod(m);

Java

import java.util.Scanner;
public class Main {
public static void main(String[] args) {
Scanner scanner = new Scanner(System.in);
int n = scanner.nextInt();
int m = scanner.nextInt();
int a = 1;
int b = 0;

Python

n, m = map(int, input().split())
a, b = 1, 0
print(0, end=" ")
for i in range(2, n + 1):
c = (i - 1) * (a + b) % m
print(c, end=" ")
a, b = b, c
print()

Stars and Bars

Resources
cp-algo
Medium

Well documented article.

Stars and Bars is a useful method in combinatorics that involves grouping indistinguishable objects into distinguishable boxes. The number of ways to put nn indistinguishable objects into kk distinguishable boxes is:

(n+k1n)=(n+k1k1)\binom{n+k-1}{n}=\binom{n+k-1}{k-1}

The second binomial coefficient from above can be derived from the property of binomial coefficients: (nk)=(nnk)\binom{n}{k}=\binom{n}{n-k}.

Let's take a look at a particular example for n=3n=3 and k=2k=2 that has 44 possibilities. As the name implies, the visualization is usually done with stars separated into groups by bars:

||\bigstar \bigstar \bigstar \\ |\bigstar |\bigstar \bigstar \\ |\bigstar \bigstar| \bigstar \\ \bigstar \bigstar \bigstar || \\

As you've probably noticed, there can be empty boxes - when we put all the stars in the first or second box. There may be cases in which the all the boxes should be non-empty. In that case, the number of ways to put nn indistinguishable objects into kk distinguishable non-empty boxes is: (n1k1)\binom{n-1}{k-1}

Focus Problem – try your best to solve this problem before continuing!

Explanation

For this problem we should think the other way around: let's say that the kk colors from which we choose are in fact boxes and instead of choosing nn marbles we put them in the respective boxes. The problem has the restriction that we should pick at least one marble of all kinds, which means in our new perspective that all the boxes should be non-empty. Thus, the answer is obtained by the second formula: (n1k1)\binom{n-1}{k-1}

Implementation

Time Complexity: O(TK)\mathcal{O}(T \cdot K)

C++

#include <iostream>
/** @return n choose k, computed naively since the problem has no overflow */
long long comb(int n, int k) {
if (k > n - k) { k = n - k; }
long long ret = 1;
for (int i = 0; i < k; i++) {
// this is done instead of *= for divisibility issues
ret = ret * (n - i) / (i + 1);
}

Python

from math import comb
ans = []
for _ in range(int(input())):
marble_num, color_num = [int(i) for i in input().split()]
ans.append(str(comb(marble_num - 1, color_num - 1)))
print("\n".join(ans))

Expected Value

Resources
Brilliant

lots of pure math examples of expected value + explanations

Brilliant

lots of pure math examples of linearity of expectation + explanations

CF

a good blog about the expected value

An expected value is the theoretical mean of a probability distribution. A random variable is used to represent a possible probability distribution. Let XX be a random variable and P(X=x)P(X = x) be the probability that the result of the random variable XX is xx. Then, the expected value of XX, denoted as E[X]E[X], is

xxP(X=x)\sum_x x \cdot P(X = x)

For example, let XX be the probability distribution of a fair 6-sided die. P(X=x)P(X = x) is 16\frac{1}{6} for 1x61 \leq x \leq 6. Using the formula, we get E[X]=216=72E[X] = \frac{21}{6} = \frac{7}{2}.

Focus Problem – try your best to solve this problem before continuing!

Explanation

Let XX represent the probability distribution of the maximum number of candies a child gets. To get E[X]E[X], we need P(X=x)P(X = x) for 1xk1 \leq x \leq k. (xk)n(\frac{x}{k})^n is the probability that each child gets at most xx candies. To get P(X=x)P(X = x), we must subtract out the probability that each child gets strictly less than xx candies, which is (x1k)n(\frac{x - 1}{k})^n.

Therefore, P(X=x)=(xk)n(x1k)nP(X = x) = (\frac{x}{k})^n - (\frac{x - 1}{k})^n, allowing us to calculate E[X]E[X].

Implementation

Time Complexity: O(nk)\mathcal{O}(nk) using naive exponentiation or O(klogn)\mathcal{O}(k \log n) using binary exponentiation.

C++

#include <cmath>
#include <iomanip>
#include <iostream>
using std::cout;
using std::endl;
int main() {
int n;
int k;

Python

n, k = [int(i) for i in input().split()]
expected_max = 0
for c in range(1, k + 1):
expected_max += c * ((c / k) ** n - ((c - 1) / k) ** n)
print(f"{expected_max:.6f}")

Linearity of Expectation

Linearity of expectation states that E[X+Y]=E[X]+E[Y]E[X + Y] = E[X] + E[Y] no matter if XX and YY are independent of each other. For example, if on a certain day Alice goes to the gym with a probability of 110\frac{1}{10} and her husband Bob goes to the gym with a probability of 310\frac{3}{10}, the expected number of visits to the gym on a certain day amongst the couple is 110+310=25\frac{1}{10} + \frac{3}{10} = \frac{2}{5}. This works even though the decisions of one person may affect the other.

This can be generalized to a sequence of random variables X1,X2,,XnX_1, X_2, \dots, X_n and arbitrary constants c1,c2,,cnc_1, c_2, \dots, c_n:

E[i=1nciXi]=i=1nciE[Xi]E\left[\sum_{i = 1}^{n} c_iX_i \right] = \sum_{i = 1}^{n}c_i \cdot E[ X_i]

Focus Problem – try your best to solve this problem before continuing!

Explanation

While this may seem impossible at first glance given the amount of different sequences of operations, linearity of expectation allows us to break down the expected value into the sum of smaller parts.

We can break this down by decomposing the initial random variable into the sums of various indicator random variables. Let XuX_u be a variable that indicates whether node uu was explicitly marked for removal. Since it's basically a boolean, it can only take on either 00 or 11, which also means that the probability that it's 11 is equal to its expected value.

Now, let's consider how an operation affects node uu.

  1. Either it chooses a node that can't reach uu, and nothing happens as far as uu is concerned.
  2. It chooses uu itself, making the indicator variable 11.
  3. It chooses another node that can reach uu, making the indicator variable 00.

Since all nodes have an equal chance of being chosen, the chance that Xu=1X_u=1 is 1au\frac{1}{a_u}, where aua_u is the number of nodes that can reach uu including uu itself.

We put the expected values of all the indicator variables back together, giving us our final answer.

Implementation

Time Complexity: O(N3)\mathcal{O}(N^3) using a naive DFS or BFS.

C++

#include <iomanip>
#include <iostream>
#include <vector>
using std::cout;
using std::endl;
using std::vector;
int main() {
int n;

Python

n = int(input())
adj = [[] for _ in range(n)]
for i in range(n):
for v, c in enumerate(input()):
if c == "1":
adj[i].append(v)
# reach_from[i] = the # of nodes that have a path to i
reach_from = [0 for _ in range(n)]

Expected Products

Linearity of expectation deals with E[X+Y]E[X + Y], but what about E[XY]E[X \cdot Y]?

E[XY]=E[X]E[Y]E[X \cdot Y] = E[X] \cdot E[Y] if XX and YY are independent from each other. We can reconsider the example of a fair 6-sided die to show that E[X2]E[X]2E[X^2] \neq E[X]^2. We know that E[X]=72E[X] = \frac{7}{2}, so E[X]2=7272=494E[X]^2 = \frac{7}{2} \cdot \frac{7}{2} = \frac{49}{4}.

On the other hand,

E[X2]=xx2P(X=x)=1+4+9+16+25+366=916.E[X^2] = \sum_x x^2 \cdot P(X = x) = \frac{1 + 4 + 9 + 16 + 25 + 36}{6} = \frac{91}{6}.

Problems

StatusSourceProblem NameDifficultyTags
CSESEasy
Show TagsCombinatorics
CSESEasy
Show TagsCombinatorics
CFEasy
Show TagsCombinatorics
CFEasy
Show TagsBinary Search, Combinatorics
BronzeEasy
Show TagsCombinatorics
Bubble CupNormal
Show TagsCombinatorics
CSESNormal
Show TagsBitwise, Combinatorics
CFNormal
Show TagsCombinatorics
ACNormal
Show TagsCombinatorics
GoldNormal
Show TagsCombinatorics
GoldNormal
Show TagsBitset, PIE
CFNormal
Show TagsCombinatorics, DP
GoldNormal
Show TagsCombinatorics, Prefix Sums
CFHard
Show TagsCombinatorics, DP
GoldHard
Show TagsBinary Search, Combinatorics, Math, Probability

Module Progress:

Join the USACO Forum!

Stuck on a problem, or don't understand a module? Join the USACO Forum and get help from other competitive programmers!

PrevNext