Central Limit Theorem - Rule of thumb for repeated sampling

Question

My question was inspired by this post which concerns some of the myths and misunderstandings surrounding the Central Limit Theorem. I was asked a question by a colleague once and I couldn't offer an adequate response/solution.

My colleague's question: Statisticians often cleave to rules of thumb for the sample size of each draw (e.g., $n = 30$, $n = 50$, $n = 100$, etc.) from a population. But is there a rule of thumb for the number of times we must repeat this process?

I replied that if we were to repeat this process of taking random draws of "30 or more" (rough guideline) from a population say "thousands and thousands" of times (iterations), then the histogram of sample means will tend towards something Gaussian-like. To be clear, my confusion is not related to the number of measurements drawn, but rather the number of times (iterations) required to achieve normality. I often describe this as some theoretical process we repeat ad infinitum.

Below this question is a quick simulation in R. I sampled from the exponential distribution. The first column of the matrix X holds the 10,000 sample means, with each mean having a sample size of 2. The second column holds another 10,000 sample means, with each mean having a sample size of 4. This process repeats for columns 3 and 4 for $n = 30$ and $n = 100$, respectively. I then produced for histograms. Note, the only thing changing between the plots is the sample size, not the number of times we calculate the sample mean. Each calculation of the sample mean for a given sample size is repeated 10,000 times. We could, however, repeat this procedure 100,000 times, or even 1,000,000 times.

Questions:

(1) Is there any criteria for the number of repetitions (iterations) we must conduct to observe normality? I could try 1,000 iterations at each sample size and achieve a reasonably similar result.

(2) Is it tenable for me to conclude that this process is assumed to be repeated thousands or even millions of times? I was taught that the number of times (repetitions/iterations) is not relevant. But maybe there was a rule of thumb before the gift of modern computing power. Any thoughts?

pop <- rexp(100000, 1/10)               # The mean of the exponential distribution is 1/lambda
X <- matrix(ncol = 4, nrow = 10000)     # 10,000 repetitions

samp_sizes <- c(2, 4, 30, 100)

for (j in 1:ncol(X)) {
  for (i in 1:nrow(X)) {
    X[i, j] <- mean(sample(pop, size = samp_sizes[j]))
  }
}

par(mfrow = c(2, 2))

for (j in 1:ncol(X)) {
  hist(X[ ,j], 
       breaks = 30, 
       xlim = c(0, 30), 
       col = "blue", 
       xlab = "", 
       main = paste("Sample Size =", samp_sizes[j]))
}

Answer 1

To facilitate accurate discussion of this issue, I am going to give a mathematical account of what you are doing. Suppose you have an infinite matrix $\mathbf{X} \equiv [X_{i,j} | i \in \mathbb{Z}, j \in \mathbb{Z} ]$ composed of IID random variables from some distribution with mean $\mu$ and finite variance $\sigma^2$ that is not a normal distribution:$^\dagger$

$$X_{i,j} \sim \text{IID Dist}(\mu, \sigma^2)$$

In your analysis you are forming repeated independent iterations of sample means based on a fixed sample size. If you use a sample size of $n$ and take $M$ iterations then you are forming the statistics $\bar{X}_n^{(1)},...,\bar{X}_n^{(M)}$ given by:

$$\bar{X}_n^{(m)} \equiv \frac{1}{n} \sum_{i=1}^n X_{i,m} \quad \quad \quad \text{for } m = 1,...,M.$$

In your output you show histograms of the outcomes $\bar{X}_n^{(1)},...,\bar{X}_n^{(M)}$ for different values of $n$. It is clear that as $n$ gets bigger, we get closer to the normal distribution.

Now, in terms of "convergence to the normal distribution" there are two issues here. The central limit theorem says that the true distribution of the sample mean will converge towards the normal distribution as $n \rightarrow \infty$ (when appropriately standardised). The law of large numbers says that your histograms will converge towards the true underlying distribution of the sample mean as $M \rightarrow \infty$. So, in those histograms we have two sources of "error" relative to a perfect normal distribution. For smaller $n$ the true distribution of the sample mean is further away from the normal distribution, and for smaller $M$ the histogram is further away from the true distribution (i.e., contains more random error).

---

How big does $n$ need to be? The various "rules of thumb" for the requisite size of $n$ are not particularly useful in my view. It is true that some textbooks propagate the notion that $n=30$ is sufficient to ensure that the sample mean is well approximated by the normal distribution. The truth is that the "required sample size" for good approximation by the normal distribution is not a fixed quantity --- it depends on two factors: the degree to which the underlying distribution departs from the normal distribution; and the required level of accuracy needed for the approximation.

The only real way to determine the appropriate sample size required for an "accurate" approximation by the normal distribution is to have a look at the convergence for a range of underlying distributions. The kinds of simulations you are doing are a good way to get a sense of this.

---

How big does $M$ need to be? There are some useful mathematical results showing the rate of convergence of an empirical distribution to the true underlying distribution for IID data. To give a brief account of this, let us suppose that $F_n$ is the true distribution function for the sample mean with $n$ values, and define the empirical distribution of the simulated sample means as:

$$\hat{F}_n (x) \equiv \frac{1}{M} \sum_{m=1}^M \mathbb{I}(\bar{X}_n^{(m)} \leqslant x) \quad \quad \quad \text{for } x \in \mathbb{R}.$$

It is trivial to show that $M \hat{F}_n(x) \sim \text{Bin}(M, F_n(x))$, so the "error" between the true distribution and the empirical distribution at any point $x \in \mathbb{R}$ has zero mean, and has variance:

$$\mathbb{V} (\hat{F}_n(x) - F_n(x)) = \frac{F_n(x) (1-F_n(x))}{M}.$$

It is fairly simple to use standard confidence interval results for the binomial distribution to get an appropriate confidence intervale for the error in the simulated estimation of the distribution of the sample mean.

---

$^\dagger$ Of course, it is possible to use a normal distribution, but that is not very interesting because convergence to normality is already achieved with a sample size of one.

Answer 2

I think it may be helpful to think about your question a bit differently. Suppose that $X\sim F_X$ where $F_X$ is any arbitrary distribution, and let $\sigma^2 = Var(X)$. Now suppose I draw iid $X_1,\dots,X_n \sim F_X$, and let $\bar{X}_n = \frac{1}{n}\sum X_i$.

The CLT says that under very weak assumptions, $\bar{X}_n \xrightarrow{d} N(\mu,\sigma^2/n)$ as $n$ gets arbitrarily large. Now suppose that for a fixed $n$, I observe $\bar{X}_{n1},\dots,\bar{X}_{nK}$ where for each $k$, I sample iid $X_{1k},\dots,X_{nk} \sim F_X$ and build $\bar{X}_{nk}$. But this is the exact same as sampling $\bar{X}_{ni}$ from the distribution $F_{\bar{X}_n}$. Your question can thus posed as follows:

What is the distribution $F_{\bar{X}_n}$, and in particular, is it normal?

The answer is no, and I'll focus on your exponential example. We can understand this problem by literally considering the sampling distribution of $\bar{X}_n$ given iid $X_1,\dots,X_n \sim Exp(\gamma)$. Note that $Exp(\gamma) = \text{Gamma}(\alpha=1,\gamma)$, and so $\sum X_i \sim \text{Gamma}(n,\gamma)$ and thus

$$\frac{1}{n}\sum X_i \sim \text{Gamma}(n,\gamma/n)$$

As it turns out, for $n$ reasonably large, this distribution is very similar to a Normal distribution, but it will never be a normal distribution for any finite $n$ (the above is exactly what distribution it is!). What you did by replicating was simply drawing from this distribution and plotting (indeed, try plotting these and youll get the same result!). Depending on the distribution of $X_i$, the distribution of $\bar{X}_n$ can be anything.

What the CLT says is that as $n$ goes to infinity, $\bar{X}_n$ will converge to a normal distribution, and similarly, $\text{Gamma}(n,\gamma/n)$ (or any $F_{\bar{X}_n}$ where $X$ satisfies the requisite requirements for CLT to kick in) will asymptotically equal a normal distribution.

EDIT

In response to your comments, maybe there's a misunderstand somewhere. It's helpful to emphasize that we can think of $\bar{X}_n$ as a random variable itself (often we think of it as the mean and thus a constant, but this is not true!). The point is that the random variable $\bar{X}_n$ that is the sample mean of $X_1,\dots,X_n \sim F_X$, and the random variable $Y \sim F_{\bar{X}_n}$ are the exact same random variable. So by drawing $K$ iid draws of $X_1,\dots,X_n \sim F_X$ and calculating $\bar{X}_n$, you're doing the equivalent of $K$ draws from $F_{\bar{X}_n}$. At the end of the day, regardless of whether $K = 100,1000,100000,\dots$, youre just drawing $K$ times from $F_{\bar{X}_n}$. So what is your goal here? Are you asking at what point does the empirical cdf of $K$ draws accurately represent the cdf of $F_{\bar{X}_N}$? Well forget about anything about sample means in that case, and simply ask how many times do I need to draw some random variable $W \sim F$ such that the empirical cdf $\hat{F}_n$ is 'approximately' $F$. Well there's a whole literature on that, and two basic results are (see the wiki link on empirical cdfs for more):

1. By the Glivenko-Cantelli theorem, $\hat{F}_n$ uniformly converges to $F$ almost surely.

2. By Donsker's theorem, The empirical process $\sqrt{n}(\hat{F}_n -F)$ converges in distribution to a mean-zero Gaussian process.

What you are doing with your histograms in your post is really estimating the density (not the CDF) given $K$ draws. Histograms are a (discrete) example of kernel density estimation (KDE). There's a similar literature on KDEs, and again, you have properties like the sample KDE will converge to the true underlying density as you gather more draws (ie $K\to\infty$). It should be noted that histograms don't converge to the true density unless you also let the bin width go to zero, and this is one reason why kernel approaches are preferred: they allow smoothness and similar properties. But at the end of the day, what you can say is the following:

For a fixed $n$, drawing iid $X_1,\dots,X_n$ and considering the random variable $\frac{1}{n}\sum_{X_i}$ is equivalent to considering the random variable with distribution $F_{\bar{X}_n}$. For any $K$ draws from $F_{\bar{X}_n}$, you can estimate the CDF (empirical CDF) and/or estimate the density (two approach are histogram or KDE). In either case, as $K\to\infty$, these two esimates will converge to the true CDF/density of the random variable $\bar{X}_n$, but these will never be the normal CDF/desntiy for any fixed $n$. However, as you let $n\to\infty$, $\bar{X}_n$ is asymptotically normal (under suitable conditions), and similarly, the CDF/density will also become normal. If you take $n\to\infty$, and then $K\to\infty$, then you will get the cdf/density of a normal rv.