marimo-team/marimo-learn
3
1# /// script2# requires-python = ">=3.10"3# dependencies = [4# "marimo",5# "matplotlib==3.10.8",6# "scipy==1.17.1",7# "numpy==2.4.3",8# "plotly==5.18.0",9# "wigglystuff==0.2.37",10# ]11# ///12 13import marimo14 15__generated_with = "0.18.4"16app = marimo.App(width="medium", app_title="Central Limit Theorem")17 18 19@app.cell(hide_code=True)20def _(mo):21 mo.md(r"""22 # Central Limit Theorem23 24 _This notebook is a computational companion to ["Probability for Computer Scientists"](https://chrispiech.github.io/probabilityForComputerScientists/en/part4/clt/), by Stanford professor Chris Piech._25 26 The central limit theorem is honestly mind-blowing — it's like magic that no matter what distribution you start with, the sampling distribution of means approaches a normal distribution as sample size increases.27 28 Mathematically, if we have:29 30 $X_1, X_2, \ldots, X_n$ as independent, identically distributed random variables with:31 32 - Mean: $\mu$33 - Variance: $\sigma^2 < \infty$34 35 Then as $n \to \infty$:36 37 $$\sqrt{n}\left(\frac{1}{n}\sum_{i=1}^{n}X_i - \mu\right) \xrightarrow{d} \mathcal{N}(0, \sigma^2)$$38 39 > _Note:_ The above LaTeX derivation is included as a reference. Credit for this formulation goes to the original source linked at the top of the notebook.40 """)41 return42 43 44@app.cell(hide_code=True)45def _(mo):46 mo.md(r"""47 ## Central Limit Theorem Statement48 49 There are two ways to state the central limit theorem:50 51 ### Sum Version52 53 Let $X_1, X_2, \dots, X_n$ be independent and identically distributed random variables. The sum of these random variables approaches a normal distribution as $n \rightarrow \infty$:54 55 $\sum_{i=1}^{n}X_i \sim \mathcal{N}(n \cdot \mu, n \cdot \sigma^2)$56 57 Where $\mu = E[X_i]$ and $\sigma^2 = \text{Var}(X_i)$. Since each $X_i$ is identically distributed, they share the same expectation and variance.58 59 ### Average Version60 61 Let $X_1, X_2, \dots, X_n$ be independent and identically distributed random variables. The average of these random variables approaches a normal distribution as $n \rightarrow \infty$:62 63 $\frac{1}{n}\sum_{i=1}^{n}X_i \sim \mathcal{N}(\mu, \frac{\sigma^2}{n})$64 65 Where $\mu = E[X_i]$ and $\sigma^2 = \text{Var}(X_i)$.66 67 The CLT is incredible because it applies to almost any distribution (as long as it has a finite mean and variance), regardless of its shape.68 """)69 return70 71 72@app.cell(hide_code=True)73def _(mo):74 mo.md(r"""75 ## Central Limit Theorem Intuition76 77 Let's explore what happens when you add random variables together. For example, what if we add 100 different uniform random variables?78 79 ```python80 from random import random81 82 def add_100_uniforms():83 total = 084 for i in range(100):85 # returns a sample from uniform(0, 1)86 x_i = random()87 total += x_i88 return total89 ```90 91 The value returned by this function will be a random variable. Click the button below to run the function and observe the resulting value of total:92 """)93 return94 95 96@app.cell(hide_code=True)97def _(mo):98 run_button = mo.ui.run_button(label="Run add_100_uniforms()")99 100 run_button.center()101 return (run_button,)102 103 104@app.cell(hide_code=True)105def _(mo, random, run_button):106 def add_100_uniforms():107 total = 0108 for i in range(100):109 # returns a sample from uniform(0, 1)110 x_i = random.random() 111 total += x_i112 return total113 114 # Display the result when the button is clicked115 if run_button.value:116 uniform_result = add_100_uniforms()117 display = mo.md(f"**total**: {uniform_result:.5f}")118 else:119 display = mo.md("")120 121 display122 return (add_100_uniforms,)123 124 125@app.cell(hide_code=True)126def _(mo):127 mo.md(r"""128 What does total look like as a distribution? Let's calculate total many times and visualize the histogram of values it produces.129 """)130 return131 132 133@app.cell(hide_code=True)134def _(mo):135 # Simulation control136 run_simulation_button = mo.ui.button(137 value=0, 138 on_click=lambda value: value + 1, 139 label="Run 10,000 more samples", 140 kind="warn"141 )142 143 run_simulation_button.center()144 return (run_simulation_button,)145 146 147@app.cell(hide_code=True)148def _(add_100_uniforms, go, mo, np, run_simulation_button, stats, time):149 # store the results150 def get_simulation_results():151 if not hasattr(get_simulation_results, "results"):152 get_simulation_results.results = []153 get_simulation_results.last_button_value = -1 # track button clicks154 return get_simulation_results155 156 # grab the results157 sim_storage = get_simulation_results()158 simulation_results = sim_storage.results159 160 # Check if button was clicked (value changed)161 if run_simulation_button.value != sim_storage.last_button_value:162 # Update the last seen button value163 sim_storage.last_button_value = run_simulation_button.value164 165 with mo.status.spinner(title="Running simulation...") as progress_status:166 sim_count = 10000167 new_results = []168 for _ in mo.status.progress_bar(range(sim_count)):169 sim_result = add_100_uniforms()170 new_results.append(sim_result)171 time.sleep(0.0001) # tiny pause172 173 simulation_results.extend(new_results)174 175 progress_status.update(f"✅ Added {sim_count:,} samples (total: {len(simulation_results):,})")176 177 if simulation_results:178 # Numbers179 mean = np.mean(simulation_results)180 std_dev = np.std(simulation_results)181 182 theoretical_mean = 100 * 0.5 # = 50183 theoretical_variance = 100 * (1/12) # = 8.33...184 theoretical_std = np.sqrt(theoretical_variance) # ≈ 2.89185 186 # should be 10k times the click number (mainly for the y-axis label)187 total_samples = run_simulation_button.value * 10000188 189 fig = go.Figure()190 191 # histogram of samples192 fig.add_trace(go.Histogram(193 x=simulation_results,194 histnorm='probability density',195 name='Sum Distribution',196 marker_color='royalblue',197 opacity=0.7198 ))199 200 x_vals = np.linspace(min(simulation_results), max(simulation_results), 1000)201 y_vals = stats.norm.pdf(x_vals, theoretical_mean, theoretical_std)202 203 fig.add_trace(go.Scatter(204 x=x_vals,205 y=y_vals,206 mode='lines',207 name='Normal approximation',208 line=dict(color='red', width=2)209 ))210 211 fig.add_vline(212 x=mean, 213 line_dash="dash", 214 line_width=1.5,215 line_color="green",216 annotation_text=f"Sample Mean: {mean:.2f}",217 annotation_position="top right"218 )219 220 # some notes221 fig.add_annotation(222 x=0.02, y=0.95,223 xref="paper", yref="paper",224 text=f"Sum of 100 Uniform(0,1) variables<br>" +225 f"Sample size: {total_samples:,}<br>" +226 f"Sample mean: {mean:.2f} (expected: {theoretical_mean})<br>" +227 f"Sample std: {std_dev:.2f} (expected: {theoretical_std:.2f})<br>" +228 f"According to CLT: Normal({theoretical_mean}, {theoretical_variance:.2f})",229 showarrow=False,230 align="left",231 bgcolor="white",232 opacity=0.8233 )234 235 fig.update_layout(236 title=f'Distribution of Sum of 100 Uniforms (Click #{run_simulation_button.value})',237 xaxis_title='Values',238 yaxis_title=f'Probability Density ({total_samples:,} runs)',239 template='plotly_white',240 height=500241 )242 243 # show244 histogram = mo.ui.plotly(fig)245 else:246 histogram = mo.md("Click the button to run the simulation!")247 248 # display249 histogram250 return251 252 253@app.cell(hide_code=True)254def _(mo):255 mo.md(r"""256 That is interesting! The sum of 100 independent uniforms looks normal. Is that a special property of uniforms? No! It turns out to work for almost any type of distribution (as long as the distribution has finite mean and variance).257 258 - Sum of 40 $X_i$ where $X_i \sim \text{Beta}(a = 5, b = 4)$? Normal.259 - Sum of 90 $X_i$ where $X_i \sim \text{Poisson}(\lambda = 4)$? Normal.260 - Sum of 50 dice-rolls? Normal.261 - Average of 10000 $X_i$ where $X_i \sim \text{Exp}(\lambda = 8)$? Normal.262 263 For any distribution, the sum or average of a sufficiently large number of independent, identically distributed random variables will be approximately normally distributed.264 """)265 return266 267 268@app.cell(hide_code=True)269def _(mo):270 mo.md(r"""271 ## Continuity Correction272 273 When using the Central Limit Theorem with discrete random variables (like a Binomial or Poisson), we need to apply a continuity correction. This is because we're approximating a discrete distribution with a continuous one (normal).274 275 The continuity correction involves adjusting the boundaries in probability calculations by ±0.5 to account for the discrete nature of the original variable.276 277 You should use a continuity correction any time your normal is approximating a discrete random variable. The rules for a general continuity correction are the same as the rules for the [binomial-approximation continuity correction](http://marimo.app/https://github.com/marimo-team/learn/blob/main/probability/14_binomial_distribution.py).278 279 In our example above, where we added 100 uniforms, a continuity correction isn't needed because the sum of uniforms is continuous. However, in examples with dice or other discrete distributions, a continuity correction would be necessary.280 """)281 return282 283 284@app.cell(hide_code=True)285def _(mo):286 mo.md(r"""287 ## Examples288 289 Let's work through some practical examples to see how the Central Limit Theorem is applied.290 """)291 return292 293 294@app.cell(hide_code=True)295def _(mo):296 mo.md(r"""297 ### Example 1: Dice Game298 299 > _Note:_ The following application demonstrates the practical use of the Central Limit Theorem. The mathematical derivation is based on concepts from ["Probability for Computer Scientists"](https://chrispiech.github.io/probabilityForComputerScientists/en/part2/clt/) by Chris Piech.300 301 Let's solve a fun probability problem: You roll a 6-sided die 10 times and let $X$ represent the total value of all 10 dice: $X = X_1 + X_2 + \dots + X_{10}$. You win if $X \leq 25$ or $X \geq 45$. What's your probability of winning?302 303 For a single die roll $X_i$, we know:304 - $E[X_i] = 3.5$305 - $\text{Var}(X_i) = \frac{35}{12}$306 307 **Solution Approach:**308 309 This is where the Central Limit Theorem shines! Since we're summing 10 independent, identically distributed random variables, we can approximate this sum with a normal distribution $Y$:310 311 $Y \sim \mathcal{N}(10 \cdot E[X_i], 10 \cdot \text{Var}(X_i)) = \mathcal{N}(35, 29.2)$312 313 Now calculating our winning probability:314 315 $P(X \leq 25 \text{ or } X \geq 45) = P(X \leq 25) + P(X \geq 45)$316 317 Since we're approximating a discrete distribution with a continuous one, we apply a continuity correction:318 319 $\approx P(Y < 25.5) + P(Y > 44.5) = P(Y < 25.5) + [1 - P(Y < 44.5)]$320 321 Converting to standard normal form:322 323 $\approx \Phi\left(\frac{25.5 - 35}{\sqrt{29.2}}\right) + \left[1 - \Phi\left(\frac{44.5 - 35}{\sqrt{29.2}}\right)\right]$324 325 $\approx \Phi(-1.76) + [1 - \Phi(1.76)]$326 327 $\approx 0.039 + (1 - 0.961) \approx 0.078$328 329 So your chance of winning is about 7.8% — not great odds, but that's probability for you!330 """)331 return332 333 334@app.cell(hide_code=True)335def _(create_dice_game_visualization, fig_to_image, mo):336 # Display visualization337 dice_game_fig = create_dice_game_visualization()338 dice_game_image = mo.image(fig_to_image(dice_game_fig), width="100%")339 340 dice_explanation = mo.md(341 r"""342 **Understanding the Visualization:**343 344 This graph shows our dice game in action. The blue bars represent the exact probability distribution for summing 10 dice, while the red curve shows our normal approximation from the Central Limit Theorem.345 346 I've highlighted the winning regions in orange:347 - The left region where $X \leq 25$348 - The right region where $X \geq 45$349 350 Together these regions cover about 7.8% of the total probability.351 352 What's fascinating here is how closely the normal curve approximates the actual discrete distribution — this is the Central Limit Theorem working its magic, even with just 10 random variables.353 """354 )355 356 mo.vstack([dice_game_image, dice_explanation])357 return358 359 360@app.cell(hide_code=True)361def _(mo):362 mo.md(r"""363 ### Example 2: Algorithm Runtime Estimation364 365 > _Note:_ The following derivation demonstrates the practical application of the Central Limit Theorem for experimental design. The mathematical approach is based on concepts from ["Probability for Computer Scientists"](https://chrispiech.github.io/probabilityForComputerScientists/en/part2/clt/) by Chris Piech.366 367 Here's a practical problem I encounter in performance testing: You've developed a new algorithm and want to measure its average runtime. You know the variance is $\sigma^2 = 4 \text{ sec}^2$, but need to estimate the true mean runtime $t$.368 369 The question: How many test runs do you need to be 95% confident your estimated mean is within ±0.5 seconds of the true value?370 371 Let $X_i$ represent the runtime of the $i$-th test (for $1 \leq i \leq n$).372 373 **Solution:**374 375 We need to find $n$ such that:376 377 $0.95 = P\left(-0.5 \leq \frac{\sum_{i=1}^n X_i}{n} - t \leq 0.5\right)$378 379 The Central Limit Theorem tells us that as $n$ increases, the sample mean approaches a normal distribution. Let's standardize this to work with the standard normal distribution:380 381 $Z = \frac{\left(\sum_{i=1}^n X_i\right) - n\mu}{\sigma \sqrt{n}} = \frac{\left(\sum_{i=1}^n X_i\right) - nt}{2 \sqrt{n}}$382 383 Rewriting our probability constraint in terms of $Z$:384 385 $0.95 = P\left(-0.5 \leq \frac{\sum_{i=1}^n X_i}{n} - t \leq 0.5\right) = P\left(\frac{-0.5 \sqrt{n}}{2} \leq Z \leq \frac{0.5 \sqrt{n}}{2}\right)$386 387 Using the properties of the standard normal CDF:388 389 $0.95 = \Phi\left(\frac{\sqrt{n}}{4}\right) - \Phi\left(-\frac{\sqrt{n}}{4}\right) = 2\Phi\left(\frac{\sqrt{n}}{4}\right) - 1$390 391 Solving for $\Phi\left(\frac{\sqrt{n}}{4}\right)$:392 393 $0.975 = \Phi\left(\frac{\sqrt{n}}{4}\right)$394 395 Using the inverse CDF:396 397 $\Phi^{-1}(0.975) = \frac{\sqrt{n}}{4}$398 399 $1.96 = \frac{\sqrt{n}}{4}$400 401 $n = 61.4$402 403 Rounding up, we need 62 test runs to achieve our desired confidence interval — a practical result we can immediately apply to our testing protocol.404 """)405 return406 407 408@app.cell(hide_code=True)409def _(create_algorithm_runtime_visualization, fig_to_image, mo):410 # Display visualization411 runtime_fig = create_algorithm_runtime_visualization()412 runtime_image = mo.image(fig_to_image(runtime_fig), width="100%")413 414 runtime_explanation = mo.md(415 r"""416 **Visualization Explanation:**417 418 The graph illustrates how the standard error of the mean (SEM) decreases as the number of trials increases. The standard error is calculated as $\frac{\sigma}{\sqrt{n}}$.419 420 - When we conduct 62 trials, the standard error is approximately 0.254 seconds.421 - With a 95% confidence level, this gives us a margin of error of about ±0.5 seconds (1.96 × 0.254 ≈ 0.5).422 - The shaded region shows how the confidence interval narrows as the number of trials increases.423 424 This demonstrates why 62 trials are sufficient to meet our requirements of estimating the mean runtime within ±0.5 seconds with 95% confidence.425 """426 )427 428 mo.vstack([runtime_image, runtime_explanation])429 return430 431 432@app.cell(hide_code=True)433def _(mo):434 mo.md(r"""435 ## Interactive CLT Explorer436 437 Let's explore how the Central Limit Theorem works with different underlying distributions. You can select a distribution type and see how the distribution of the sample mean changes as the sample size increases.438 """)439 return440 441 442@app.cell(hide_code=True)443def _(controls):444 controls445 return446 447 448@app.cell(hide_code=True)449def _(450 distribution_type,451 fig_to_image,452 mo,453 np,454 plt,455 run_explorer_button,456 sample_size,457 sim_count_slider,458 stats,459):460 # Run simulation when button is clicked461 if run_explorer_button.value:462 # Set distribution parameters based on selection463 if distribution_type.value == "uniform":464 dist_name = "Uniform(0, 1)"465 # For uniform(0,1): mean = 0.5, variance = 1/12466 true_mean = 0.5467 true_var = 1/12468 469 # generate samples470 def generate_sample():471 return np.random.uniform(0, 1, sample_size.value)472 473 elif distribution_type.value == "exponential":474 rate = 1.0475 dist_name = f"Exponential(λ={rate})"476 # For exponential(λ): mean = 1/λ, variance = 1/λ²477 true_mean = 1/rate478 true_var = 1/(rate**2)479 480 def generate_sample():481 return np.random.exponential(1/rate, sample_size.value)482 483 elif distribution_type.value == "binomial":484 n_param, p = 10, 0.3485 dist_name = f"Binomial(n={n_param}, p={p})"486 # For binomial(n,p): mean = np, variance = np(1-p)487 true_mean = n_param * p488 true_var = n_param * p * (1-p)489 490 def generate_sample():491 return np.random.binomial(n_param, p, sample_size.value)492 493 elif distribution_type.value == "poisson":494 rate = 3.0495 dist_name = f"Poisson(λ={rate})"496 # For poisson(λ): mean = λ, variance = λ497 true_mean = rate498 true_var = rate499 500 def generate_sample():501 return np.random.poisson(rate, sample_size.value)502 503 # Generate the simulation data using a spinner for progress504 with mo.status.spinner(title="Running simulation...") as explorer_progress:505 sample_means = []506 original_samples = []507 508 # Run simulations509 for _ in mo.status.progress_bar(range(sim_count_slider.value)):510 sample = generate_sample()511 512 # Store the first simulation's individual values for visualizing original distribution513 if len(original_samples) < 1000: # limit to prevent memory issues514 original_samples.extend(sample)515 516 # sample mean517 sample_means.append(np.mean(sample))518 519 # progress520 explorer_progress.update(f"✅ Completed {sim_count_slider.value:,} simulations")521 522 # Create visualization523 explorer_fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 5))524 525 # Original distribution histogram526 ax1.hist(original_samples, bins=30, density=True, alpha=0.7, color='royalblue')527 ax1.set_title(f"Original Distribution: {dist_name}")528 529 # Theoretical mean line530 ax1.axvline(x=true_mean, color='red', linestyle='--', 531 label=f'True Mean = {true_mean:.3f}')532 533 ax1.set_xlabel("Value")534 ax1.set_ylabel("Density")535 ax1.legend()536 537 # Sample means histogram and normal approximation538 sample_mean_mean = np.mean(sample_means)539 sample_mean_std = np.std(sample_means)540 expected_std = np.sqrt(true_var / sample_size.value) # CLT prediction541 542 ax2.hist(sample_means, bins=30, density=True, alpha=0.7, color='forestgreen',543 label=f'Sample Size = {sample_size.value}')544 545 # Normal approximation from CLT546 explorer_x = np.linspace(min(sample_means), max(sample_means), 1000)547 explorer_y = stats.norm.pdf(explorer_x, true_mean, expected_std)548 ax2.plot(explorer_x, explorer_y, 'r-', linewidth=2, label='CLT Normal Approximation')549 550 # Add mean line551 ax2.axvline(x=true_mean, color='purple', linestyle='--',552 label=f'True Mean = {true_mean:.3f}')553 554 ax2.set_title(f"Distribution of Sample Means\n(CLT Prediction: N({true_mean:.3f}, {true_var/sample_size.value:.5f}))")555 ax2.set_xlabel("Sample Mean")556 ax2.set_ylabel("Density")557 ax2.legend()558 559 # Add CLT description560 explorer_fig.text(0.5, 0.01, 561 f"Central Limit Theorem: As sample size increases, the distribution of sample means approaches\n" +562 f"a normal distribution with mean = {true_mean:.3f} and variance = {true_var:.3f}/{sample_size.value} = {true_var/sample_size.value:.5f}",563 ha='center', fontsize=10, bbox=dict(facecolor='white', alpha=0.8))564 565 plt.tight_layout(rect=[0, 0.05, 1, 1])566 567 # Display plot568 explorer_image = mo.image(fig_to_image(explorer_fig), width="100%")569 else:570 explorer_image = mo.md("Click the 'Run Simulation' button to see how the Central Limit Theorem works.")571 572 explorer_image573 return574 575 576@app.cell(hide_code=True)577def _(mo):578 mo.md(r"""579 ## 🤔 Test Your Understanding580 581 /// details | What is the shape of the distribution of the sum of many independent random variables?582 The sum of many independent random variables approaches a normal distribution, regardless of the shape of the original distributions (as long as they have finite mean and variance). This is the essence of the Central Limit Theorem.583 ///584 585 /// details | If $X_1, X_2, \dots, X_{100}$ are IID random variables with $E[X_i] = 5$ and $Var(X_i) = 9$, what is the distribution of their sum?586 By the Central Limit Theorem, the sum $S = X_1 + X_2 + \dots + X_{100}$ follows a normal distribution with:587 588 - Mean: $E[S] = 100 \cdot E[X_i] = 100 \cdot 5 = 500$589 - Variance: $Var(S) = 100 \cdot Var(X_i) = 100 \cdot 9 = 900$590 591 Therefore, $S \sim \mathcal{N}(500, 900)$, or equivalently $S \sim \mathcal{N}(500, 30^2)$.592 ///593 594 /// details | When do you need to apply a continuity correction when using the Central Limit Theorem?595 You need to apply a continuity correction when you're using the normal approximation (through CLT) for a discrete random variable.596 597 For example, when approximating a binomial or Poisson distribution with a normal distribution, you should adjust boundaries by ±0.5 to account for the discrete nature of the original variable. This makes the approximation more accurate.598 ///599 600 /// details | If $X_1, X_2, \dots, X_{n}$ are IID random variables, how does the variance of their sample mean $\bar{X} = \frac{1}{n}\sum_{i=1}^{n}X_i$ change as $n$ increases?601 The variance of the sample mean decreases as the sample size $n$ increases. Specifically:602 603 $Var(\bar{X}) = \frac{Var(X_i)}{n}$604 605 This means that as we take more samples, the sample mean becomes more concentrated around the true mean of the distribution. This is why larger samples give more precise estimates.606 ///607 608 /// details | Why is the Central Limit Theorem so important in statistics?609 The Central Limit Theorem is foundational in statistics because:610 611 1. It allows us to make inferences about population parameters using sample statistics, regardless of the population's distribution.612 2. It explains why the normal distribution appears so frequently in natural phenomena.613 3. It enables the construction of confidence intervals and hypothesis tests for means, even when the underlying population distribution is unknown.614 4. It justifies many statistical methods that assume normality, even when working with non-normal data, provided the sample size is large enough.615 616 In essence, the CLT provides the theoretical justification for much of statistical inference.617 ///618 """)619 return620 621 622@app.cell(hide_code=True)623def _(mo):624 mo.md(r"""625 ## Appendix (helper code and functions)626 """)627 return628 629 630@app.cell631def _():632 import marimo as mo633 return (mo,)634 635 636@app.cell(hide_code=True)637def _():638 from wigglystuff import TangleSlider639 return640 641 642@app.cell(hide_code=True)643def _():644 # Import libraries645 import numpy as np646 import matplotlib.pyplot as plt647 from scipy import stats648 import io649 import base64650 import random651 import time652 import plotly.graph_objects as go653 import plotly.io as pio654 return base64, go, io, np, plt, random, stats, time655 656 657@app.cell(hide_code=True)658def _(base64, io):659 from matplotlib.figure import Figure660 661 # Helper function to convert matplotlib figures to images662 def fig_to_image(fig):663 buf = io.BytesIO()664 fig.savefig(buf, format='png', bbox_inches='tight')665 buf.seek(0)666 img_str = base64.b64encode(buf.getvalue()).decode('utf-8')667 return f"data:image/png;base64,{img_str}"668 return (fig_to_image,)669 670 671@app.cell(hide_code=True)672def _(np, plt, stats):673 def create_dice_game_visualization():674 """Create a visualization for the dice game example."""675 # Parameters676 n_dice = 10677 dice_values = np.arange(1, 7) # 1 to 6678 679 # Theoretical values680 single_die_mean = np.mean(dice_values) # 3.5681 single_die_var = np.var(dice_values) # 35/12682 683 # Sum distribution parameters684 sum_mean = n_dice * single_die_mean685 sum_var = n_dice * single_die_var686 sum_std = np.sqrt(sum_var)687 688 # Possible outcomes for the sum of 10 dice689 min_sum = n_dice * min(dice_values) # 10690 max_sum = n_dice * max(dice_values) # 60691 sum_values = np.arange(min_sum, max_sum + 1)692 693 # Create figure694 fig, ax = plt.subplots(figsize=(10, 6))695 696 # Calculate PMF through convolution697 # For one die698 single_pmf = np.ones(6) / 6699 700 sum_pmf = single_pmf.copy()701 for _ in range(n_dice - 1):702 sum_pmf = np.convolve(sum_pmf, single_pmf)703 704 # Plot the PMF705 ax.bar(sum_values, sum_pmf, alpha=0.7, color='royalblue', label='Exact PMF')706 707 # Normal approximation708 x = np.linspace(min_sum - 5, max_sum + 5, 1000)709 y = stats.norm.pdf(x, sum_mean, sum_std)710 ax.plot(x, y, 'r-', linewidth=2, label='Normal Approximation')711 712 # Win conditions (x ≤ 25 or x ≥ 45)713 win_region_left = sum_values <= 25714 win_region_right = sum_values >= 45715 716 # Shade win regions717 ax.bar(sum_values[win_region_left], sum_pmf[win_region_left], 718 color='darkorange', alpha=0.7, label='Win Region')719 ax.bar(sum_values[win_region_right], sum_pmf[win_region_right], 720 color='darkorange', alpha=0.7)721 722 # Calculate win probability723 win_prob = np.sum(sum_pmf[win_region_left]) + np.sum(sum_pmf[win_region_right])724 725 # Add vertical lines for critical values726 ax.axvline(x=25.5, color='red', linestyle='--', linewidth=1.5, label='Critical Points')727 ax.axvline(x=44.5, color='red', linestyle='--', linewidth=1.5)728 729 # Add mean line730 ax.axvline(x=sum_mean, color='green', linestyle='--', linewidth=1.5, 731 label=f'Mean = {sum_mean}')732 733 # Text box with relevant information734 textstr = '\n'.join((735 f'Number of dice: {n_dice}',736 f'Sum Mean: {sum_mean}',737 f'Sum Std Dev: {sum_std:.2f}',738 f'Win Probability: {win_prob:.4f}',739 f'CLT Approximation: {0.078:.4f}'740 ))741 props = dict(boxstyle='round', facecolor='wheat', alpha=0.5)742 ax.text(0.05, 0.95, textstr, transform=ax.transAxes, fontsize=10,743 verticalalignment='top', bbox=props)744 745 # Formatting746 ax.set_xlabel('Sum of 10 Dice')747 ax.set_ylabel('Probability')748 ax.set_title('Central Limit Theorem: Dice Game Example')749 ax.legend()750 ax.grid(alpha=0.3)751 752 plt.tight_layout()753 plt.gca()754 return fig755 return (create_dice_game_visualization,)756 757 758@app.cell(hide_code=True)759def _(np, plt):760 def create_algorithm_runtime_visualization():761 """Create a visualization for the algorithm runtime example."""762 # Parameters763 variance = 4 # σ² = 4 sec²764 std_dev = np.sqrt(variance) # σ = 2 sec765 confidence_level = 0.95766 z_score = 1.96 # for 95% confidence767 target_error = 0.5 # ±0.5 seconds768 769 # Calculate n needed for desired precision770 n_required = int(np.ceil((z_score * std_dev / target_error) ** 2)) # ≈ 62771 772 n_values = np.arange(1, 100)773 774 # standard error775 standard_errors = std_dev / np.sqrt(n_values)776 777 # margin of error778 margins_of_error = z_score * standard_errors779 780 # Create figure781 fig, ax = plt.subplots(figsize=(10, 6))782 783 # standard error vs sample size plot784 ax.plot(n_values, standard_errors, 'b-', linewidth=2, label='Standard Error of Mean')785 786 # Plot margin of error vs sample size787 ax.plot(n_values, margins_of_error, 'r--', linewidth=2, 788 label=f'{confidence_level*100}% Margin of Error')789 790 ax.axvline(x=n_required, color='green', linestyle='-', linewidth=1.5,791 label=f'Required n = {n_required}')792 793 ax.axhline(y=target_error, color='purple', linestyle='--', linewidth=1.5,794 label=f'Target Error = ±{target_error} sec')795 796 # Shade the region below target error797 ax.fill_between(n_values, 0, target_error, alpha=0.2, color='green')798 799 # intersection point800 ax.plot(n_required, target_error, 'ro', markersize=8)801 ax.annotate(f'({n_required}, {target_error} sec)',802 xy=(n_required, target_error),803 xytext=(n_required + 5, target_error + 0.1),804 arrowprops=dict(facecolor='black', shrink=0.05, width=1))805 806 # Text box with appropriate information807 textstr = '\n'.join((808 f'Algorithm Variance: {variance} sec²',809 f'Standard Deviation: {std_dev} sec',810 f'Confidence Level: {confidence_level*100}%',811 f'Z-score: {z_score}',812 f'Target Error: ±{target_error} sec',813 f'Required Sample Size: {n_required}'814 ))815 props = dict(boxstyle='round', facecolor='wheat', alpha=0.5)816 ax.text(0.05, 0.95, textstr, transform=ax.transAxes, fontsize=10,817 verticalalignment='top', bbox=props)818 819 # Formatting820 ax.set_xlabel('Sample Size (n)')821 ax.set_ylabel('Error (seconds)')822 ax.set_title('Sample Size Determination for Algorithm Runtime Estimation')823 ax.set_xlim(0, 100)824 ax.set_ylim(0, 2)825 ax.legend()826 ax.grid(alpha=0.3)827 828 plt.tight_layout()829 return fig830 return (create_algorithm_runtime_visualization,)831 832 833@app.cell(hide_code=True)834def _(mo):835 mo.md(r"""836 ## Summary837 838 The Central Limit Theorem is truly one of the most remarkable ideas in all of statistics. It tells us that when we add up many independent random variables, their sum will follow a normal distribution, regardless of what the original distributions looked like. This is why we see normal distributions so often in real life – many natural phenomena are the result of numerous small, independent factors adding up.839 840 What makes the CLT so powerful is its universality. Whether we're working with dice rolls, measurement errors, or stock market returns, as long as we have enough independent samples, their average or sum will be approximately normal. For sums, the distribution will be $\mathcal{N}(n\mu, n\sigma^2)$, and for averages, it's $\mathcal{N}(\mu, \frac{\sigma^2}{n})$.841 842 The CLT gives us the foundation for confidence intervals, hypothesis testing, and many other statistical tools. Without it, we'd have a much harder time making sense of data when we don't know the underlying population distribution. Just remember that if you're working with discrete distributions, you'll need to apply a continuity correction to get more accurate results.843 844 Next time you see a normal distribution in data, think about the Central Limit Theorem – it might be the reason behind that familiar bell curve!845 """)846 return847 848 849@app.cell(hide_code=True)850def _(mo):851 # controls for the interactive explorer852 distribution_type = mo.ui.dropdown(853 options=["uniform", "exponential", "binomial", "poisson"],854 value="uniform",855 label="Distribution Type"856 )857 858 sample_size = mo.ui.slider(859 start =1,860 stop =100,861 step=1,862 value=30,863 label="Sample Size (n)"864 )865 866 sim_count_slider = mo.ui.slider(867 start =100,868 stop =10000,869 step=100,870 value=1000,871 label="Number of Simulations"872 )873 874 run_explorer_button = mo.ui.run_button(label="Run Simulation", kind="warn")875 876 controls = mo.hstack([877 mo.vstack([distribution_type, sample_size, sim_count_slider]),878 run_explorer_button879 ], justify='space-around')880 return (881 controls,882 distribution_type,883 run_explorer_button,884 sample_size,885 sim_count_slider,886 )887 888 889if __name__ == "__main__":890 app.run()891