Markdown

title: “Sampling Distributions and Central Limit Theorem Analysis” output: html_document —

Step 1: Population Distribution & Parameters

First, we load the necessary packages and dataset. We treat the NCbirths dataset as our population and visualize the distribution of birth weights (BirthWeightGm).

# Load required libraries
library(Stat2Data)
library(dplyr)
library(ggplot2)

# Load dataset
data("NCbirths")

# Plot histogram of the population variable
ggplot(NCbirths, aes(x = BirthWeightGm)) +
  geom_histogram(binwidth = 250, fill = "skyblue", color = "black") +
  labs(
    title = "Population Distribution of Birth Weights",
    x = "Birth Weight (grams)",
    y = "Count"
  ) +
  theme_minimal()

# Calculate and store population parameters
mu <- mean(NCbirths$BirthWeightGm, na.rm = TRUE)
sigma <- sd(NCbirths$BirthWeightGm, na.rm = TRUE)

mu
## [1] 3295.619
sigma
## [1] 632.977

Interpretation: The population distribution of birth weights is unimodal and slightly left-skewed. The population mean (\(\mu\)) is approximately 3295.62 grams, and the population standard deviation (\(\sigma\)) is approximately 632.98 grams.


Step 2: Simulating the Sampling Distribution

Next, we draw 1,000 random samples of size \(n = 30\) with replacement from the population and plot the sampling distribution of the sample means.

# Set sample size and random seed
n <- 30
set.seed(123)

# Draw 1,000 samples of size n=30
sample_means <- replicate(1000, {
  sample_data <- sample(NCbirths$BirthWeightGm[!is.na(NCbirths$BirthWeightGm)], size = n, replace = TRUE)
  mean(sample_data)
})

# Data frame for plotting
sim_data <- data.frame(sample_means = sample_means)

# Plot histogram of sample means
ggplot(sim_data, aes(x = sample_means)) +
  geom_histogram(binwidth = 25, fill = "lightgreen", color = "black") +
  labs(
    title = "Sampling Distribution of Sample Means (n = 30)",
    x = "Sample Mean Birth Weight (grams)",
    y = "Count"
  ) +
  theme_minimal()

# Summary statistics of simulated sample means
mean_sample_means <- mean(sample_means)
sd_sample_means <- sd(sample_means)

mean_sample_means
## [1] 3302.995
sd_sample_means
## [1] 114.0927

Interpretation: The sampling distribution of sample means is unimodal and symmetric, forming a bell curve. Compared to the population distribution, the center remains essentially the same (\(\text{mean} \approx\) 3303 grams vs. \(\mu \approx\) 3295.62 grams), but the spread is significantly narrower (\(\text{SD} \approx\) 114.09 grams vs. \(\sigma \approx\) 632.98 grams).


Step 3: Theoretical CLT Model vs. Simulation

According to the Central Limit Theorem (CLT), the standard error is \(\text{SE} = \frac{\sigma}{\sqrt{n}}\). We compare the theoretical normal model with our simulation.

# Compute Standard Error
SE <- sigma / sqrt(n)
SE
## [1] 115.5653
# Density histogram with superimposed CLT normal curve
ggplot(sim_data, aes(x = sample_means)) +
  geom_histogram(aes(y = after_stat(density)), binwidth = 25, fill = "lightgreen", color = "black") +
  stat_function(
    fun = dnorm,
    args = list(mean = mu, sd = SE),
    color = "red",
    linewidth = 1
  ) +
  labs(
    title = "Simulated Sample Means with CLT Normal Overlay",
    x = "Sample Mean Birth Weight (grams)",
    y = "Density"
  ) +
  theme_minimal()

Interpretation: The theoretical CLT model matches the simulated sampling distribution extremely well. The theoretical mean (\(\mu =\) 3295.62) is virtually identical to the simulated mean (3303), and the theoretical standard error (\(\text{SE} =\) 115.57) closely matches the simulated standard deviation (114.09). The red normal curve overlays the green histogram seamlessly.


Step 4: Finding Probability with pnorm()

We compute the theoretical probability that a random sample of \(n = 30\) babies has a mean birth weight greater than \(3,600\text{ grams}\) using pnorm(), and compare it to the proportion from our simulation.

# Theoretical probability using CLT
prob_clt <- pnorm(3600, mean = mu, sd = SE, lower.tail = FALSE)
prob_clt
## [1] 0.004221207
# Empirical proportion from simulation
prop_sim <- mean(sample_means > 3600)
prop_sim
## [1] 0.003

Interpretation: The theoretical probability from the CLT model is approximately 0.0042, and the empirical proportion from the simulation is 0.003. These values are nearly identical because 1,000 random samples provide a very close approximation to the theoretical normal curve, with minor differences due to random sampling variability in the simulation.


Step 5: Finding Cutoffs with qnorm()

We find the cutoff value above which \(90\%\) of sample means fall (the 10th percentile) using qnorm(), and compare it with the empirical 10th percentile from our simulation.

# Theoretical cutoff using CLT (10th percentile)
cutoff_clt <- qnorm(0.10, mean = mu, sd = SE)
cutoff_clt
## [1] 3147.516
# Empirical 10th percentile from simulation
cutoff_sim <- quantile(sample_means, 0.10)
cutoff_sim
##      10% 
## 3150.536

Interpretation: The theoretical 10th percentile from qnorm() is 3147.5 grams, and the empirical 10th percentile from the simulation is 3150.5 grams. In context, \(90\%\) of random samples of 30 babies will have a mean birth weight above approximately \(3,334.1\text{ grams}\).