Home

Understanding Complement and Substitute Dynamics: A Rewrite of Existing Research

Author's Photo

This article bridges the gap between theoretical research and practical application by exploring and rewriting the methodology of the paper Extracting Complements and Substitutes from Sales Data: A Network Perspective (Tian et al., 2021). The goal is to make the math and intuition accessible to applied data scientists, and to show how the ideas translate directly into a large-scale Spark implementation. Sincere gratitude to the authors for their foundational contributions.

Two core concepts Complements are products that are bought together — hot dogs and hot dog buns, shampoo and conditioner. Demand for one creates demand for the other.

Substitutes are products that can replace each other — Brand A tomatoes and Brand B tomatoes, Pepsi and Coca-Cola. A customer buys one instead of the other.

1   Why This Is Useful

Knowing which products complement or substitute each other unlocks a wide range of business decisions that would otherwise rely on intuition or expensive field studies.

Application How complement / substitute knowledge helps Example
Shelf & aisle placement Place complements near each other to trigger additional purchases. Salad tomatoes next to cucumbers → both sell more.
Product bundling Bundle strong complements at a slight discount to increase basket size. Hot dogs + hot dog buns bundle deal.
Inventory management When a product is out of stock, surface its substitutes automatically. Brand A pasta out of stock → recommend Brand B pasta.
Targeted promotions Discount a complement to drive sales of a high-margin partner product. Discount coffee filters to boost premium coffee sales.
Personalised recommendations "Customers also bought…" powered by complement scores. Add pasta to cart → suggest pasta sauce and parmesan.
Assortment optimisation Identify redundant substitutes and rationalise the product range. Five near-identical own-brand cereals → keep the top two.
Pricing strategy Raise the price of a substitute only when its pair is on promotion. Raise Brand B price when Brand A is on sale.
Trend detection Track how complement/substitute relationships shift over time. Oat milk becoming a substitute for dairy milk over 12 months.

2   The Network Approach

Rather than relying on price elasticity data (which requires controlled experiments) or word embeddings (which are hard to interpret), this paper models sales data directly as a bipartite network — a graph with two distinct types of nodes.

What is a bipartite network? A bipartite network has two sets of nodes where edges only connect nodes from different sets — never within the same set. Here:
  • Transaction nodes — one node per shopping trip / basket.
  • Product nodes — one node per unique product (TCIN / SKU).
  • Edges — a transaction node is connected to a product node if that product was purchased in that transaction.
Bipartite product-purchase network
A bipartite product-purchase network. Blue squares = transactions; red circles = products. An edge means "this product was in this basket."

The Biadjacency Matrix

The entire network is encoded in a single matrix called the biadjacency matrix \(A^{(b)}\):

\[ A^{(b)} = (A_{li}) \;\in\; \{0,1\}^{n_t \times n_p} \]

Every row is one basket; every column is one product. A 1 in cell \((l, i)\) means "product \(i\) was in basket \(l\)." This sparse binary matrix is the only input the entire method needs.

Why a network instead of a table? A flat co-occurrence count table tells you how often two products appear together, but not whether that frequency is surprising given how popular each product is individually. A network model lets us apply statistical null models to answer: "Is this co-occurrence more (or less) than we would expect by chance?"

3   Key Assumptions

The method rests on six assumptions that translate economic intuition into measurable network properties. The first four define what complements and substitutes look like in the data; the last two control for noise.

Assumptions about product relationships

  1. Complements co-occur more than expected. If two products are genuine complements, they will appear in the same basket significantly more often than a random model would predict — even after accounting for how popular each product is on its own.
  2. Stronger complements co-occur more frequently. The degree of complementarity is positively correlated with the relative frequency of co-purchase. Hot dogs and hot dog buns (bought together almost every time) are stronger complements than, say, milk and bread.
  3. Substitutes share the same complements but rarely co-occur. Brand A tomatoes and Brand B tomatoes are both complements of cucumber and lettuce, but customers rarely buy both brands in the same basket — they pick one.
  4. Stronger substitutes have more similar complement profiles. The more alike two products' sets of complements are, the more interchangeable those products are. This is the key insight that lets us measure substitutability without any price data.

Assumptions about noise

  1. Noise does not change the community structure. Random co-purchases (e.g., two popular products that happen to appear in the same large basket) may create spurious local connections, but they will not distort the overall groups (communities) of complements and substitutes. This justifies using community detection as a validation step.
  2. Noise can be modelled and removed. The null models described in the next section capture the expected co-occurrence due to popularity and basket-size effects alone. Any co-occurrence beyond what the null model predicts is treated as a genuine signal.
Why these assumptions matter for the implementation Assumption 3 is what makes the substitute detection two-stage: first find complements, then find products that share complements but avoid each other. This is exactly the logic in the Spark notebook — pre_supplementary_pairs is filtered to same-department pairs with low co-occurrence, then further filtered to those sharing at least one common complement.

4   Null Models — Deciding What Is Significant

A null model is a randomised version of the network that preserves certain structural properties (like how popular each product is) but removes any genuine product relationships. By comparing the real network to the null model, we can identify co-occurrences that are statistically surprising.

The key quantity: common neighbours \(cn_{ij}\) For any two products \(i\) and \(j\), their common neighbours \(cn_{ij}\) is simply the number of transactions in which both products appear together: \[ cn_{ij} = \sum_{l=1}^{n_t} A_{li} \cdot A_{lj} \] If \(cn_{ij}\) is much higher than the null model predicts → likely complements.
If \(cn_{ij}\) is much lower than the null model predicts → likely substitutes.

Two null models are proposed, in increasing order of sophistication:

4.1   Bipartite Erdős–Rényi (ER) Model

The ER model is the simpler of the two. It asks: "If each product independently appeared in each transaction with some fixed probability, how many shared transactions would we expect for a pair of products?"

Core idea

Each product \(i\) is assigned a probability \(p_i\) of appearing in any given transaction. This probability is estimated directly from the data:

\[ \hat{p}_i = \frac{d_i^{(p)}}{n_t} \]

Intuitively: if a product appears in 10% of all baskets, its probability is 0.10. A very popular product (milk) has a high \(p_i\); a niche product has a low \(p_i\).

Expected co-occurrence

Under this model, the probability that a single transaction contains both products \(i\) and \(j\) is simply \(p_i \cdot p_j\) (since they appear independently). The number of shared transactions \(cn_{ij}\) therefore follows a Binomial distribution:

\[ X_{ij} \;\sim\; \text{Binomial}(n_t,\; p_i p_j) \]

For large \(n_t\), the Central Limit Theorem lets us approximate this with a Normal distribution:

\[ \mu_{ij} = n_t p_i p_j \qquad \sigma_{ij}^2 = n_t p_i p_j (1 - p_i p_j) \]

Significance thresholds

We declare a co-occurrence significantly more than expected (→ complement candidate) if:

\[ cn_{ij} \;>\; \mu_{ij} + \Phi^{-1}(1 - \alpha_m)\,\sigma_{ij} \]

And significantly less than expected (→ substitute candidate) if:

\[ cn_{ij} \;<\; \mu_{ij} - \Phi^{-1}(1 - \alpha_l)\,\sigma_{ij} \]

Here \(\Phi^{-1}(\cdot)\) is the inverse CDF of the standard Normal distribution, \(\alpha_m\) is the significance level for "more" (e.g., 0.01), and \(\alpha_l\) is the significance level for "less" (e.g., 0.20).

Limitation of the ER model The ER model assumes every transaction has the same size. In reality, some customers buy 2 items and others buy 30. A large basket inflates \(cn_{ij}\) for all pairs of products in it, even if those products have nothing to do with each other. The BiCM below fixes this.

4.2   Bipartite Configuration Model (BiCM)

The BiCM is a more realistic null model. It preserves both the degree of each product node (how popular each product is) and the degree of each transaction node (how large each basket is). This means the null model accounts for the fact that a large basket will naturally contain more product pairs.

How it works conceptually

Imagine each product node has \(d_i^{(p)}\) "stubs" (half-edges) and each transaction node has \(d_l^{(t)}\) stubs. The BiCM randomly connects these stubs, but only allows edges between transaction nodes and product nodes (never product–product or transaction–transaction). The resulting random network has the same degree sequence as the real data, but no genuine product relationships.

Probability of a shared transaction

Under the BiCM, the probability that products \(i\) and \(j\) both appear in transaction \(l\) is:

\[ p_{ilj} = \frac{d_i^{(p)} \cdot d_l^{(t)} \cdot d_j^{(p)} \cdot (d_l^{(t)} - 1)}{m^2} \]

Notice that \(p_{ilj}\) varies per transaction: a large basket (high \(d_l^{(t)}\)) contributes more to the expected co-occurrence than a small basket. This is the key improvement over the ER model.

Expected co-occurrence \(\mu_{ij}\) — the key formula

Summing \(p_{ilj}\) over all transactions gives the expected number of shared transactions:

\[ \mu_{ij} = \sum_{l=1}^{n_t} p_{ilj} = \frac{d_i^{(p)} \cdot d_j^{(p)}}{m} \cdot \frac{\langle d^{(t)2} \rangle - \langle d^{(t)} \rangle}{\langle d^{(t)} \rangle} \]
Why this formula is a computational breakthrough Look at the right-hand side carefully. It splits into two parts:
  • Product-pair part: \(\dfrac{d_i^{(p)} \cdot d_j^{(p)}}{m}\) — depends only on the two products being compared.
  • Transaction moment part: \(\dfrac{\langle d^{(t)2} \rangle - \langle d^{(t)} \rangle}{\langle d^{(t)} \rangle}\) — a single scalar constant computed once from all transactions, the same for every product pair.
This means we can compute \(\mu_{ij}\) for billions of product pairs in Spark with a simple join and multiplication — no per-pair loop over transactions needed. In the notebook this constant is called moment_simplified.

Poisson approximation

The number of shared transactions \(cn_{ij}\) is a sum of independent Bernoulli variables (one per transaction), each with probability \(p_{ilj}\). This is called a Poisson binomial distribution. For sparse networks where each \(p_{ilj}\) is small, Le Cam's theorem guarantees it can be well approximated by a simpler Poisson distribution:

\[ cn_{ij} \;\approx\; Y_{ij} \;\sim\; \text{Poisson}(\mu_{ij}) \]

This approximation is valid because real retail networks are sparse — most products appear in only a small fraction of all transactions, so \(p_{ilj}\) is tiny for most pairs.

Significance testing with the Poisson CDF

Let \(F_{ij}(y)\) be the cumulative distribution function of \(\text{Poisson}(\mu_{ij})\):

\[ F_{ij}(y) = e^{-\mu_{ij}} \sum_{k=0}^{y} \frac{\mu_{ij}^k}{k!} \]

We then apply two one-sided tests:

Test Condition Interpretation Notebook threshold
Significantly more \(1 - F_{ij}(cn_{ij}) < \alpha_m\) Co-occurrence is surprisingly high → complement candidate p_value > 0.85
Significantly less \(F_{ij}(cn_{ij}) < \alpha_l\) Co-occurrence is surprisingly low → substitute candidate p_value < 0.15
Important: minimum co-occurrence filter The Poisson CDF test alone is not sufficient at retail scale. A pair that co-occurred just once with an expected mean of 0.01 gets \(F(1, 0.01) \approx 0.9999\), trivially passing the complement threshold. At 837K products and 2B co-occurring pairs, this inflates the complement count to over 1 billion. A minimum co-occurrence filter (e.g., count ≥ 5) must be applied before the significance test to ensure only meaningfully co-occurring pairs are classified as complements.

5   Building the Complement & Substitute Networks

The significance tests produce two binary (unweighted) unipartite networks — graphs where every node is a product and edges indicate a statistically significant relationship.

Step 1 — Two intermediate networks from the significance tests

MatrixMeaningEntry = 1 when…
\(A^{(m)}\) "More" network Products \(i\) and \(j\) co-occur significantly more than expected
\(A^{(l)}\) "Less" network Products \(i\) and \(j\) co-occur significantly less than expected

Step 2 — Derive the complement and substitute networks

Complement network \(A^{(c)}\)

By Assumption 1, complements are simply the pairs that co-occur significantly more than expected. So the complement network is identical to \(A^{(m)}\):

\[ A^{(c)} = A^{(m)} \]

Substitute network \(A^{(s)}\)

By Assumption 3, substitutes must satisfy two conditions simultaneously:

  1. They co-occur less than expected — captured by \(A^{(l)}\).
  2. They share at least one common complement — captured by checking whether \((A^{(m)})^T A^{(m)}\) has a non-zero entry for that pair.

Both conditions are combined with an element-wise (Hadamard) product:

\[ A^{(s)} = \mathbf{I}_{\{(A^{(m)})^T A^{(m)} > 0\}} \;\odot\; A^{(l)} \]

The indicator matrix \(\mathbf{I}_{\{(A^{(m)})^T A^{(m)} > 0\}}\) acts as a gate: it is 1 for a pair \((i,j)\) only if there exists at least one product \(k\) such that both \(i\) and \(j\) are complements of \(k\). Without this gate, any two unpopular products that happen to rarely co-occur would be labelled substitutes — which is clearly wrong.

Concrete example Brand A tomatoes and Brand B tomatoes:
  • Both are complements of cucumber, lettuce, and olive oil → \((A^{(m)})^T A^{(m)} > 0\) ✓
  • Customers rarely buy both brands in the same basket → \(A^{(l)} = 1\) ✓
  • Therefore \(A^{(s)}_{ij} = 1\) → they are substitutes ✓

At this stage both networks are unweighted — they only tell us whether a relationship exists, not how strong it is. The next two sections add weights.

6   Measuring Complementarity — simo(i, j)

Knowing that two products are complements is useful, but knowing how strongly complementary they are is even more valuable. The paper proposes a family of measures derived from a weighted cosine similarity on the bipartite network. The simplest and most practical is the original measure, simo.

Intuition: a random-walk perspective

Imagine placing a unit of "signal" on product \(i\) and letting it flow one step through the bipartite network — from product \(i\) to all transactions containing \(i\), weighted by \(1/d_l^{(t)}\) (the inverse basket size). The resulting distribution over transactions is the "fingerprint" of product \(i\).

Two products are strongly complementary if their fingerprints overlap a lot — i.e., they tend to appear in the same transactions, especially small transactions where the co-occurrence is less likely to be coincidental.

The formula

\[ simo(i, j) = \frac{ \displaystyle\sum_{l=1}^{n_t} \frac{A_{li}\, A_{lj}}{d_l^{(t)}} }{ \sqrt{ \left(\displaystyle\sum_{h=1}^{n_t} \frac{A_{hi}}{d_h^{(t)}}\right) \left(\displaystyle\sum_{h=1}^{n_t} \frac{A_{hj}}{d_h^{(t)}}\right) } } \]

Breaking it down term by term

TermWhat it computesWhy it matters
\(A_{li} \cdot A_{lj}\) 1 if both products appear in transaction \(l\), else 0 Selects only transactions where both products co-occur
\(\dfrac{1}{d_l^{(t)}}\) Inverse basket size of transaction \(l\) Discounts large baskets. If a basket has 50 items, any two products in it are likely there by coincidence. A basket of 2 items is a much stronger signal of genuine complementarity.
Numerator sum \(\sum_l A_{li} A_{lj} / d_l^{(t)}\) Total "weighted co-occurrence" of products \(i\) and \(j\)
Denominator \(\sum_h A_{hi}/d_h^{(t)}\) Total weighted occurrence of product \(i\) alone Normalises by how often each product appears individually, so that popular products don't automatically get high scores with everything.
\(\sqrt{\cdots \times \cdots}\) Geometric mean of the two individual weighted occurrences Makes the score symmetric: \(simo(i,j) = simo(j,i)\), and ensures \(simo(i,i) = 1\) (a product is perfectly complementary to itself).
Range and interpretation \(simo(i,j) \in [0, 1]\). A value near 1 means the two products almost always appear together in small, focused baskets — very strong complements. A value near 0 means they rarely share transactions, or only share large baskets where the co-occurrence is likely coincidental.

From unweighted to weighted complement network

The simo score is computed for all product pairs, but only pairs already identified as complements (i.e., \(A^{(c)}_{ij} = 1\)) are assigned a non-zero weight. This is done via an element-wise product:

\[ W^{(c)} = A^{(c)} \;\odot\; simo(i,j) \]

\(A^{(c)}\) acts as a mask: it zeroes out the score for any pair that did not pass the significance test, no matter how high their raw simo value might be. This prevents noise from inflating the complement network.

Spark implementation note

In the notebook, the numerator is computed by joining the sales data on transaction_id (to find co-occurring pairs) and summing 1/d_t per pair. The denominator is computed by summing 1/d_t per product individually, then joining back. The critical fix (Bug 1) was changing F.sum(F.col("tcin_1")) to F.sum(F.col("inverse")) in the denominator — the original code was summing the product ID number instead of the inverse transaction degree.

7   Measuring Substitutability — sims(i, j)

Once we have the weighted complement network \(W^{(c)}\), we can measure how substitutable two products are. The core idea (Assumption 4) is:

Substitutes have similar complement profiles. Each product \(i\) can be described by a vector of its complementarity scores with every other product: \(\mathbf{w}_i = (W^{(c)}_{i1}, W^{(c)}_{i2}, \ldots, W^{(c)}_{in_p})\). Two products are strong substitutes if their complement vectors point in the same direction — i.e., they are complementary to the same set of products with similar strengths.

The formula — cosine similarity over complement vectors

\[ sims(i, j) = \frac{ \displaystyle\sum_{k=1}^{n_p} W^{(c)}_{ik}\, W^{(c)}_{jk} }{ \sqrt{ \left(\displaystyle\sum_{p=1}^{n_p} \left(W^{(c)}_{ip}\right)^2\right) \left(\displaystyle\sum_{p=1}^{n_p} \left(W^{(c)}_{jp}\right)^2\right) } } \]

Breaking it down

TermWhat it computesWhy it matters
\(\sum_k W^{(c)}_{ik} W^{(c)}_{jk}\) Dot product of the two complement vectors Large when both products have high complementarity scores with the same third products \(k\). This is the numerator of cosine similarity.
\(\sum_p (W^{(c)}_{ip})^2\) Squared norm of product \(i\)'s complement vector Normalises for the overall "strength" of product \(i\)'s complement relationships, so that a product with many strong complements doesn't automatically score high with everything.
\(\sqrt{\cdots \times \cdots}\) Geometric mean of the two norms Makes the score symmetric and bounded in \([0, 1]\).
Concrete example Brand A tomatoes has high simo scores with: cucumber (0.8), lettuce (0.7), olive oil (0.5).
Brand B tomatoes has high simo scores with: cucumber (0.75), lettuce (0.65), olive oil (0.45).

Their complement vectors are nearly parallel → high cosine similarity → high sims → strong substitutes. ✓

Why the θ_c threshold matters before computing sims

The paper (§3.4.2, Appendix B) explicitly warns that low-quality complementarity scores must be removed before computing sims. If near-zero simo values are included, they add noise to the complement vectors and bias the cosine similarity. The paper calibrates a threshold \(\theta_c\) at the 0.35 quantile of non-zero simo values (approximately 0.011 in their dataset). In the notebook:

theta_c = measure_complementary.approxQuantile('simo', [0.35], 0.01)[0]
measure_complementary = measure_complementary.filter(F.col('simo') >= theta_c)

From unweighted to weighted substitute network

Just as with complements, the sims score is masked by the binary substitute network \(A^{(s)}\):

\[ W^{(s)} = A^{(s)} \;\odot\; sims(i,j) \]

Only pairs that passed the significance test (low co-occurrence + shared complement) receive a non-zero substitutability score. The resulting matrix \(W^{(s)}\) is the final output: a weighted network where edge weight = strength of substitutability.

Range and interpretation \(sims(i,j) \in [0, 1]\). A value near 1 means the two products have almost identical complement profiles — they are nearly perfect substitutes. A value near 0 means their complement profiles are very different, suggesting they serve different purposes despite rarely being bought together.

8   End-to-End Pipeline Summary

The full method goes from raw transaction data to two weighted product networks in seven steps:

1
Build biadjacency matrix \(A^{(b)}\): transactions × products
→
2
Compute expected co-occurrence \(\mu_{ij}\) via BiCM formula using moment_simplified
→
3
Poisson CDF significance test Filter: count ≥ 5, then classify as complement or substitute candidate
→
4
Build \(A^{(c)}\) and \(A^{(s)}\) Substitute network requires shared complement gate
→
5
Compute simo Weighted cosine similarity on bipartite network; apply \(\theta_c\) threshold
→
6
Compute sims Cosine similarity over thresholded complement vectors
→
7
Output \(W^{(c)}\) and \(W^{(s)}\) Weighted complement and substitute networks, ready for community detection or ranking

Key parameters (from paper Appendix B)

ParameterValueMeaning
\(\alpha_m\) 0.01 Significance level for "significantly more" co-occurrence (complement test)
\(\alpha_l\) 0.20 Significance level for "significantly less" co-occurrence (substitute test)
\(\theta_c\) 0.35 quantile of non-zero simo Threshold below which complementarity scores are zeroed out before computing sims
MIN_CO_OCCURRENCE 5 (practical addition for large-scale data) Minimum observed co-occurrence before applying the Poisson CDF test

Bugs found and fixed in the Spark implementation

BugSymptomRoot causeFix
Bug 1 (Critical) simo ≈ 1e-10, sims ≈ 1e-8 F.sum(F.col("tcin_1")) summed product IDs instead of 1/d_t Changed to F.sum(F.col("inverse"))
Bug 2 1.1B complementary pairs (too many) Pairs with count=1 and tiny expected mean trivially pass the CDF test Added count ≥ 5 filter before significance test
Bug 2b 700K substitute pairs (too few) Zero-cooccurrence pairs (true substitutes) absent from co-occurrence matrix Generate same-department zero-count candidates; compute proper \(p = e^{-\mu}\)
Bug 3 Same pair repeated 10+ times in output drop_duplicates called before tcin_3 fan-out was resolved Aggregate countDistinct(tcin_3) first, then dedup
Bug 4 sims biased by noise No \(\theta_c\) threshold applied before computing sims Apply approxQuantile('simo', [0.35]) threshold per paper Appendix B
Bug 5 False substitute signals for rare product pairs Zero-count pairs assigned blanket p_value=0 regardless of \(\mu\) Compute p_value = exp(-expected_mean) per pair; filter by p_value < 0.15