> For the complete documentation index, see [llms.txt](https://zavalab.gitbook.io/statschemebook/llms.txt). Markdown versions of documentation pages are available by appending `.md` to page URLs; this page is available as [Markdown](https://zavalab.gitbook.io/statschemebook/chapter-4-estimation-for-random-variables/07-central-limit-theorem.md).

# Central Limit Theorem

* Now turn our attention to assessing accuracy of sample approx $$\hat{\mathbb{E}}\_X^S$$.
* Because MC approx is an RV, accuracy is measured in terms of its uncertainty.
* If we want uncertainty of an RV, we need to know its underlying pdf.
* Obtaining pdf of $$\hat{\mathbb{E}}\_X^S$$ can be achieved by using a powerful result in statistics known as the central limit theorem (CLT).
* Consider a sample sequence $$X\_1,X\_2,...,X\_S$$ that is i.i.d. and that has *known* expected value $$e=\mathbb{E}\[X]$$ and variance $$v^2=\mathbb{V}\[X]$$.
* CLT states that:

  $$
  \lim\_{S\to \infty}\hat{\mathbb{E}}\_X^S\sim \mathcal{N}(e,v^2/S)
  $$
* CLT states that, *regardless* of the nature of $$X$$ (e.g., Weibull, Poisson), sample approximation $$\hat{\mathbb{E}}\_X^S$$ will *always* behave as a Gaussian RV as $$S$$ increases.
* CLT states that expected value of $$\hat{\mathbb{E}}\_X^S$$ is true value $$e$$. In other words, $$\hat{\mathbb{E}}\_X^S$$ is an unbiased estimate of $$e$$.
* CLT states that the variance $$\mathbb{V}\[\hat{\mathbb{E}}\_X^S]$$ is given by $$v^2/S$$ (std dev is $$v/\sqrt{S}$$) and that this shrinks with $$S$$. In other words, $$\hat{\mathbb{E}}\_X^S$$ becomes more accurate as $$S$$ increases.
* Since we know that $$\hat{\mathbb{E}}\_X^S\sim \mathcal{N}(e,v^2/S)$$ then we have pdf for $$\hat{\mathbb{E}}\_X^S$$ and thus we can compute quantities of interest for it.
* For example, we can compute $$1-\alpha$$ confidence intervals:

  $$
  \mathbb{E}\_X\in \left\[\hat{\mathbb{E}}*X^S\pm z*{\alpha/2}v/\sqrt{S}\right].
  $$
* Note that this interval shrinks as $$S$$ increases.
* CLT assumes that we know variance of $$X$$; as a result, we need to know this quantity to compute confidence interval.
* In practice, one often replaces this quantity with its MC approx (sample variance):

  $$
  \hat{v}^2=\frac{1}{S-1}\sum\_{\omega\in\mathcal{S}}(x\_\omega-\hat{\mathbb{E}}\_X^S)^2.
  $$

## Example: Verifying the Central Limit Theorem

{% hint style="danger" icon="flask" %}

* We draw i.i.d. samples $$X\_1,...,X\_S$$ from $$X\sim \textrm{Weibull}(1,2)$$ and compute $$\hat{\mathbb{E}}*X^S=\frac{1}{S}\sum*{\omega=1}^S X\_\omega$$.
* We repeat procedure for different sets of samples and visualize pdf of $$\bar{X}=\hat{\mathbb{E}}\_X^S$$.
  {% endhint %}

<figure><img src="https://1753452229-files.gitbook.io/~/files/v0/b/gitbook-x-prod.appspot.com/o/spaces%2FakKE106umUCSWFtTnl31%2Fuploads%2Fgit-blob-172510356759cbb490be7b663b534f7559918dbe%2Fch4_clt_weibull.png?alt=media" alt="" width="350"><figcaption><p>Pdf for 𝑋 ∼ Weibull(1, 2) (top-left) and sample average 𝑋̄ with 𝑆 = 10 (top-right), 𝑆 = 100 (bottom-left) and 𝑆 = 1000 (bottom-right).</p></figcaption></figure>

#### Python implementation

{% hint style="success" %}
Notebook [`ch4_clt_example.ipynb`](https://colab.research.google.com/github/vzavala/StatsBook/blob/master/Python/ch4_clt_example.ipynb) — open it in Colab and run the cells alongside this section.

MATLAB live script: [`ch4_clt_example.mlx`](https://github.com/vzavala/StatsBook/blob/master/Matlab-Live/ch4_clt_example.mlx)
{% endhint %}

Verifying the central limit theorem (Victor Z, UW-Madison, 2024). $$X$$ is Weibull with scale $$\beta=2$$ and shape $$\xi=1$$ (as in the slide figure).

| Slides                                       | Python                | Meaning                              |
| -------------------------------------------- | --------------------- | ------------------------------------ |
| $$S$$                                        | `S`                   | sample size of each average          |
| $$\bar{X}=\hat{\mathbb{E}}^S\_X$$            | `x_bar`               | sample averages (one per repetition) |
| $$e=\mathbb{E}\[X]$$, $$v^2=\mathbb{V}\[X]$$ | `e`, `v2`             | true mean and variance               |
| $$\mathcal{N}(e, v^2/S)$$                    | `norm(e, sqrt(v2/S))` | CLT limit                            |

`R` is the number of repetitions used to draw the histogram of $$\bar{X}$$.

**Note:** This notebook reproduces Example 4.7 (Section 4.3.3) of *Statistics for Chemical Engineers*: the distribution of the sample average of Weibull(1,2) observations.

```python
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import weibull_min, norm

xi, beta = 1, 2
X = weibull_min(xi, scale=beta)
e, v2 = X.mean(), X.var()
R = 1000

np.random.seed(0)
fig, axs = plt.subplots(2, 2, figsize=(10, 8))
x = np.linspace(0, 10, 400)
axs[0, 0].plot(x, X.pdf(x), 'black', linewidth=1.5)
axs[0, 0].set_xlabel('$x$')
axs[0, 0].set_ylabel('$f(x)$')
z = np.linspace(0, 4, 400)
for ax, S in zip((axs[0, 1], axs[1, 0], axs[1, 1]), (10, 100, 1000)):
    x_bar = X.rvs(size=(R, S)).mean(axis=1)               # R independent sample averages
    ax.hist(x_bar, bins=30, density=True, edgecolor='black', facecolor='none')
    ax.plot(z, norm.pdf(z, e, np.sqrt(v2 / S)), 'k', linewidth=1.5)  # CLT: N(e, v^2/S)
    ax.set_xlim(0, 4)
    ax.set_title(f'$S={S}$')
    ax.set_xlabel(r'$\bar{x}$')
    ax.set_ylabel(r'$f(\bar{x})$')
    print(f"S = {S:4d}:  mean(x_bar) = {x_bar.mean():.3f} (e = {e:.3f}),  "
          f"S * var(x_bar) = {S * x_bar.var():.3f} (v^2 = {v2:.3f})")
for ax in axs.ravel():
    ax.grid(True)
plt.tight_layout()
plt.show()
```

```
S =   10:  mean(x_bar) = 1.981 (e = 2.000),  S * var(x_bar) = 4.073 (v^2 = 4.000)
S =  100:  mean(x_bar) = 2.004 (e = 2.000),  S * var(x_bar) = 3.909 (v^2 = 4.000)
S = 1000:  mean(x_bar) = 2.003 (e = 2.000),  S * var(x_bar) = 3.929 (v^2 = 4.000)
```

<figure><img src="https://1753452229-files.gitbook.io/~/files/v0/b/gitbook-x-prod.appspot.com/o/spaces%2FakKE106umUCSWFtTnl31%2Fuploads%2Fgit-blob-35b219fdd3b7ef4aacd78111be57b120a177b470%2Fch4_clt_example-03-1.png?alt=media" alt="" width="520"><figcaption></figcaption></figure>

**Note:** As in the book, the sample average is Gaussian even though X is not, and its variance shrinks with S: this is the central limit theorem.