Skip to content
Computational Lab.

Probability Distributions

The formula, and what sampling actually gives you.

Beginner / probability / sampling

The exact probability mass function sits on top of a running histogram of seeded samples for the binomial, Poisson, or normal distribution. The empirical mean is checked against the theoretical one in the metrics - sampling noise shrinking in public.

Tagged: binomial, Poisson, normal, pmf

Live simulation

The accent curve is the exact pmf; the bars are what sampling gives.

Speed

The accent curve is the exact pmf; the bars are what sampling gives.

Mathematics

Three named distributions, three different stories: binomial counts successes in n trials, Poisson counts rare events per interval, normal approximates everything the CLT touches. Each has a closed-form pmf that can be drawn exactly - and then tested against sampling.

binomial: P(k) = C(n, k) p^k (1-p)^(n-k)

Poisson: P(k) = lambda^k e^{-lambda} / k!

normal: f(x) = exp(-(x-mu)^2 / 2 sigma^2) / (sigma sqrt(2 pi))

The exact pmf curve is drawn from the formulas above; the histogram accumulates seeded samples whose empirical mean is compared with the theoretical mean in the metrics. The binomial coefficients are computed in log space so n up to 20 never overflows.

How it works

Sampling per distribution: binomial sums n Bernoulli draws, Poisson uses Knuth's multiplication loop, normal uses Box-Muller. The histogram normalizes to its own maximum so the pmf overlay scales to the same axis.

Change parameters and the histogram resets - a different distribution needs a fresh count.

Implementation

src/engine/experiments/probability-distributions.ts · 180 lines · the module running above

TypeScript · fixed-timestep engine

1import { BaseExperiment } from "../core/base-experiment";
2import type { ParameterDef, ThemeColors } from "../core/types";
3
4const BINS = 40;
5
6/**
7 * Discrete distributions three ways: the exact pmf (computed from its
8 * formula) against a running histogram of seeded samples. Binomial,
9 * Poisson, and normal curves with parameters you can drag while the
10 * sampling keeps accumulating.
11 */
12export class ProbabilityDistributionsExperiment extends BaseExperiment {
13 readonly id = "probability-distributions";
14
15 private bins = new Float64Array(BINS);
16 private count = 0;
17 private domainMax = 20;
18
19 protected params(): ParameterDef[] {
20 return [
21 {
22 key: "distribution",
23 label: "Distribution",
24 options: [
25 { value: "binomial", label: "Binomial(n, p)" },
26 { value: "poisson", label: "Poisson(lambda)" },
27 { value: "normal", label: "Normal(mu, sigma)" },
28 ],
29 defaultValue: "binomial",
30 },
31 { key: "paramA", label: "n / lambda / mu", min: 1, max: 20, step: 0.5, defaultValue: 12 },
32 { key: "paramB", label: "p / - / sigma", min: 0.05, max: 1, step: 0.05, defaultValue: 0.5 },
33 { key: "rate", label: "Samples per tick", min: 10, max: 1000, step: 10, defaultValue: 100 },
34 ];
35 }
36
37 protected onReset(): void {
38 this.bins.fill(0);
39 this.count = 0;
40 }
41
42 protected onParameterChange(key: string): void {
43 if (key !== "rate") this.reset();
44 }
45
46 private domain(): number {
47 const dist = this.str("distribution");
48 const a = this.num("paramA");
49 if (dist === "binomial") return Math.min(39, Math.ceil(a) + 1);
50 if (dist === "poisson") return Math.min(39, Math.ceil(a * 2.2));
51 return 20; // normal domain: mu +- 5 sigma mapped into BINS
52 }
53
54 /** Exact probability mass at integer k (normal: bin prob over its bin). */
55 private pmf(k: number): number {
56 const dist = this.str("distribution");
57 const a = this.num("paramA");
58 const b = this.num("paramB");
59 if (dist === "binomial") {
60 const n = Math.round(a);
61 if (k > n) return 0;
62 let logChoose = 0;
63 for (let i = 1; i <= k; i++) logChoose += Math.log((n - k + i) / i);
64 return Math.exp(logChoose + k * Math.log(b) + (n - k) * Math.log(1 - b));
65 }
66 if (dist === "poisson") {
67 let logP = -a + k * Math.log(a);
68 for (let i = 1; i <= k; i++) logP -= Math.log(i);
69 return Math.exp(logP);
70 }
71 const mu = a;
72 const sigma = Math.max(0.05, b);
73 const x = k;
74 return (
75 Math.exp(-((x - mu) ** 2) / (2 * sigma * sigma)) /
76 (sigma * Math.sqrt(2 * Math.PI))
77 );
78 }
79
80 private drawSample(): number {
81 const dist = this.str("distribution");
82 const a = this.num("paramA");
83 const b = this.num("paramB");
84 if (dist === "binomial") {
85 const n = Math.round(a);
86 let hits = 0;
87 for (let i = 0; i < n; i++) if (this.rng() < b) hits += 1;
88 return hits;
89 }
90 if (dist === "poisson") {
91 // Knuth's method for moderate lambda.
92 const lambda = a;
93 const limit = Math.exp(-lambda);
94 let k = 0;
95 let p = 1;
96 do {
97 k += 1;
98 p *= this.rng();
99 } while (p > limit && k < 200);
100 return k - 1;
101 }
102 // Normal via Box-Muller, clamped into the domain.
103 let u = 0;
104 let v = 0;
105 while (u === 0) u = this.rng();
106 while (v === 0) v = this.rng();
107 const z = Math.sqrt(-2 * Math.log(u)) * Math.cos(2 * Math.PI * v);
108 return Math.min(BINS - 0.51, Math.max(-0.5, a + z * Math.max(0.05, b)));
109 }
110
111 protected onUpdate(): void {
112 this.domainMax = this.domain();
113 const rate = this.num("rate");
114 for (let i = 0; i < rate; i++) {
115 const value = this.drawSample();
116 const bin = Math.min(BINS - 1, Math.max(0, Math.round(value)));
117 this.bins[bin] += 1;
118 this.count += 1;
119 }
120 }
121
122 render(ctx: CanvasRenderingContext2D, theme: ThemeColors): void {
123 ctx.fillStyle = theme.bg;
124 ctx.fillRect(0, 0, this.width, this.height);
125 const binWidth = this.width / BINS;
126 const baseY = this.height - 20;
127 const chartH = this.height - 46;
128 let maxCount = 1;
129 for (let i = 0; i < BINS; i++) maxCount = Math.max(maxCount, this.bins[i]);
130 ctx.fillStyle = theme.fgSecondary;
131 for (let i = 0; i < BINS; i++) {
132 const h = (this.bins[i] / maxCount) * chartH;
133 ctx.fillRect(i * binWidth + 0.5, baseY - h, binWidth - 1, h);
134 }
135 // Exact pmf overlay scaled to the histogram's mass.
136 ctx.strokeStyle = theme.accent;
137 ctx.lineWidth = 2;
138 ctx.beginPath();
139 const toPixels = (prob: number) => prob * (this.count / BINS) * (chartH / maxCount) * BINS;
140 for (let i = 0; i < BINS; i++) {
141 const prob = this.pmf(i);
142 const y = baseY - Math.min(chartH, toPixels(prob));
143 if (i === 0) ctx.moveTo(i * binWidth + binWidth / 2, y);
144 else ctx.lineTo(i * binWidth + binWidth / 2, y);
145 }
146 ctx.stroke();
147 ctx.fillStyle = theme.lineStrong;
148 ctx.fillRect(0, baseY, this.width, 1);
149 ctx.font = "11px ui-monospace, monospace";
150 ctx.fillStyle = theme.fgTertiary;
151 ctx.fillText(`n = ${this.count.toLocaleString("en-US")} samples`, 10, 16);
152 }
153
154 getMetrics() {
155 // Empirical mean vs theoretical mean over the binned samples.
156 let sum = 0;
157 for (let i = 0; i < BINS; i++) sum += this.bins[i] * i;
158 const empirical = this.count > 0 ? sum / this.count : 0;
159 const dist = this.str("distribution");
160 const a = this.num("paramA");
161 const b = this.num("paramB");
162 const theoretical = dist === "binomial" ? Math.round(a) * b : dist === "poisson" ? a : a;
163 return {
164 samples: this.count.toLocaleString("en-US"),
165 distribution: dist,
166 empiricalMean: empirical.toFixed(3),
167 theoreticalMean: theoretical.toFixed(3),
168 };
169 }
170
171 describe(): string {
172 const m = this.getMetrics();
173 return `${m.samples} samples from the ${m.distribution} distribution: empirical mean ${m.empiricalMean} against the theoretical ${m.theoreticalMean}.`;
174 }
175
176 entities(): number {
177 return this.count;
178 }
179}
180

Performance

Time per tick
O(1) per sample (binomial O(n))
Space
One 40-bin histogram
Notes
Poisson with large lambda approaches the normal - switch between them at lambda = 12 and see it.

Live numbers (FPS, frame time, steps per second, microseconds per step, entity count) are in the performance overlay on the simulation above: the pulse icon in its top-right corner.