⚙️ Setup & data¶
import numpy as np, pandas as pd, matplotlib.pyplot as plt
from scipy import stats
AMBER="#d97706"; TEAL="#0d9488"; INK="#1a2138"; GRID="#e6e9f2"; PINK="#db2777"
plt.rcParams.update({"figure.facecolor":"white","axes.facecolor":"white","figure.dpi":110,"font.size":11,
"axes.edgecolor":GRID,"axes.grid":True,"grid.color":GRID,"axes.axisbelow":True,"axes.spines.top":False,
"axes.spines.right":False,"axes.titlesize":12,"axes.titleweight":"bold","legend.frameon":False})
BASE="https://raw.githubusercontent.com/johnfisher-ai/Statistics-Data-Science-AI-Visual-Book/main/data/"
try: df = pd.read_csv("../../data/ecommerce_session_logs.csv")
except FileNotFoundError: df = pd.read_csv(BASE+"ecommerce_session_logs.csv")
print("loaded:", df.shape)
df.head()
loaded: (1000, 5)
| session_id | timestamp | user_id | device_type | cart_addition | |
|---|---|---|---|---|---|
| 0 | SESS_10000 | 2026-01-01 00:38:00 | USR_22081 | Desktop | 0 |
| 1 | SESS_10001 | 2026-01-01 01:06:00 | USR_23419 | Desktop | 1 |
| 2 | SESS_10002 | 2026-01-01 00:58:00 | USR_21574 | Mobile | 0 |
| 3 | SESS_10003 | 2026-01-01 00:59:00 | USR_20249 | Mobile | 0 |
| 4 | SESS_10004 | 2026-01-01 01:42:00 | USR_20673 | Mobile | 0 |
x = df["cart_addition"]
p = x.mean()
print(f"sessions : {len(x):,}")
print(f"add-to-cart events : {int(x.sum()):,}")
print(f"baseline rate p : {p:.4f} ({p:.1%})")
print(f"Bernoulli mean = p = {p:.4f}")
print(f"Bernoulli variance = p(1-p) = {p*(1-p):.4f}")
print(f"std dev = sqrt(p(1-p)) = {np.sqrt(p*(1-p)):.4f}")
sessions : 1,000 add-to-cart events : 148 baseline rate p : 0.1480 (14.8%) Bernoulli mean = p = 0.1480 Bernoulli variance = p(1-p) = 0.1261 std dev = sqrt(p(1-p)) = 0.3551
fig,ax=plt.subplots(figsize=(5.2,3.4))
obs=[(x==0).mean(),(x==1).mean()]
ax.bar([0,1],obs,color=[GRID,AMBER],edgecolor=INK,width=0.6)
for i,v in enumerate(obs): ax.text(i,v+0.01,f"{v:.3f}",ha="center",fontweight="bold")
ax.set_xticks([0,1],["no cart (0)","add to cart (1)"]); ax.set_ylim(0,1)
ax.set_ylabel("proportion of sessions"); ax.set_title(f"Bernoulli outcomes: p = {p:.3f}")
plt.tight_layout(); plt.show()
The whole distribution is two bars: about 85% of sessions add nothing, 15% convert. That single rate p is the heartbeat of the funnel, and its variance p(1−p) is largest when p is near 0.5 and small when conversions are rare, as here.
by_dev = df.groupby("device_type")["cart_addition"].agg(["mean","count"]).sort_values("mean",ascending=False)
by_dev.columns=["conversion_rate","sessions"]
print(by_dev.round(4))
conversion_rate sessions device_type Tablet 0.2381 84 Mobile 0.1398 608 Desktop 0.1396 308
fig,ax=plt.subplots(figsize=(6,3.2))
ax.bar(by_dev.index, by_dev["conversion_rate"], color=AMBER, edgecolor=INK, width=0.6)
ax.axhline(p, color=PINK, ls="--", lw=2, label=f"overall p = {p:.3f}")
for i,(r) in enumerate(by_dev["conversion_rate"]): ax.text(i,r+0.004,f"{r:.3f}",ha="center",fontweight="bold")
ax.set_ylabel("conversion rate"); ax.set_title("Conversion rate by device"); ax.legend()
plt.tight_layout(); plt.show()
Each device is its own Bernoulli process with its own p. Wherever a segment sits clearly above or below the dashed overall rate, that is a lever: the business can chase the gap. (These are still just estimates, Chapter 59 shows how to test whether such differences are real.)
rng = np.random.default_rng(48)
# simulate a long stream of sessions; record the gap (in sessions) between conversions
stream = rng.random(2_000_000) < p
idx = np.flatnonzero(stream)
gaps = np.diff(idx) # sessions between consecutive conversions
print(f"simulated mean gap = {gaps.mean():.2f} sessions (theory 1/p = {1/p:.2f})")
print(f"P(next conversion within 5 sessions) = {1-(1-p)**5:.3f}")
print(f"P(next conversion within 10 sessions) = {1-(1-p)**10:.3f}")
simulated mean gap = 6.76 sessions (theory 1/p = 6.76) P(next conversion within 5 sessions) = 0.551 P(next conversion within 10 sessions) = 0.798
fig,ax=plt.subplots(figsize=(7,3.2))
kmax=40; ax.hist(gaps,bins=np.arange(1,kmax+1),density=True,color=AMBER,alpha=0.7,label="simulated gaps")
k=np.arange(1,kmax); ax.plot(k, stats.geom.pmf(k,p),"o-",color=TEAL,ms=4,label="Geometric(p) PMF")
ax.axvline(1/p,color=PINK,ls="--",lw=2,label=f"mean = 1/p = {1/p:.1f}")
ax.set_xlabel("sessions until next add-to-cart"); ax.set_ylabel("probability"); ax.set_title("When will the next purchase arrive?"); ax.legend()
plt.tight_layout(); plt.show()
The simulated waiting times trace the geometric PMF exactly, with an average wait of about 7 sessions between conversions. The business can now answer "how long until the next sale?" probabilistically: there is a 55% chance within 5 sessions and a 80% chance within 10. Bernoulli trials, stacked in time, became a forecast.
n_day = 200
expected = n_day*p
print(f"in a day of {n_day} sessions: expected conversions = n*p = {expected:.1f}")
print(f"std dev = sqrt(n p (1-p)) = {np.sqrt(n_day*p*(1-p)):.2f}")
sim_days = rng.binomial(n_day, p, size=100_000)
print(f"simulated daily mean = {sim_days.mean():.1f}, P(>=40 conversions) = {(sim_days>=40).mean():.3f}")
in a day of 200 sessions: expected conversions = n*p = 29.6 std dev = sqrt(n p (1-p)) = 5.02
simulated daily mean = 29.6, P(>=40 conversions) = 0.028
A single yes/no trial scales all the way up: one session is Bernoulli, a fixed batch of them is binomial, and a long run of daily totals is approximately normal (the Central Limit Theorem, Chapter 57). The same p that described one click now forecasts the whole day, with an expected 30 conversions per 200 sessions.