Lab 6: Inference for categorical data

library(tidyverse)
library(openintro)
library(infer)

Exercise 1

What are the counts within each category for the amount of days these students have texted while driving within the past 30 days?

Answer:

yrbss %>% 
  count(text_while_driving_30d) %>% 
  arrange(desc(n))

## # A tibble: 9 × 2
##   text_while_driving_30d     n
##   <chr>                  <int>
## 1 0                       4792
## 2 did not drive           4646
## 3 1-2                      925
## 4 <NA>                     918
## 5 30                       827
## 6 3-5                      493
## 7 10-19                    373
## 8 6-9                      311
## 9 20-29                    298

Exercise 2

What is the proportion of people who have texted while driving every day in the past 30 days and never wear helmets? Answer:

There are two common interpretations for this question depending on whether you are looking at the conditional proportion among non-helmet wearers (the standard method in this lab) or the joint proportion across the entire dataset. 1. Conditional Proportion (Among non-helmet wearers)

# Filter for students who never wear helmets
no_helmet <- yrbss %>% 
  filter(helmet_12m == "never") %>% 
  mutate(text_ind = ifelse(text_while_driving_30d == "30", "yes", "no"))

# Calculate proportion (excluding NAs)
no_helmet %>% 
  filter(!is.na(text_ind)) %>% 
  summarise(
    every_day = sum(text_ind == "yes"),
    total = n(),
    proportion = mean(text_ind == "yes")
  )

## # A tibble: 1 × 3
##   every_day total proportion
##       <int> <int>      <dbl>
## 1       463  6503     0.0712

yrbss %>% 
  summarise(
    prop = sum(text_while_driving_30d == "30" & helmet_12m == "never", na.rm = TRUE) / n()
  )

## # A tibble: 1 × 1
##     prop
##    <dbl>
## 1 0.0341

data('yrbss', package='openintro')
no_helmet <- yrbss %>%
  filter(helmet_12m == "never")

no_helmet <- no_helmet %>%
  mutate(text_ind = ifelse(text_while_driving_30d == "30", "yes", "no"))

Ecercise 3

What is the margin of error for the estimate of the proportion of non-helmet wearers that have texted while driving each day for the past 30 days based on this survey? Answer: The Margin of Error (\(\text{MoE}\)) for a symmetric confidence interval is half the width of the interval: \[\text{MoE} = \frac{\text{Upper Limit} - \text{Lower Limit}}{2}\]From the infer output provided:Lower Limit (\(\text{Lower CI}\)): \(0.0650\) (\(6.50\%\))Upper Limit (\(\text{Upper CI}\)): \(0.0775\) (\(7.75\%\))Calculation:\[\text{MoE} = \frac{0.0775 - 0.0650}{2} = \frac{0.0125}{2} = 0.00625 \quad (0.625\%)\]Point Estimate (\(\hat{p}\)): \(\frac{0.0775 + 0.0650}{2} = 0.07125 \quad (7.125\%)\)Result: The margin of error for this estimate is \(0.00625\) (or \(0.625\%\)).

Exercise 4

Using the infer package, calculate confidence intervals for two other categorical variables (you’ll need to decide which level to call “success”, and report the associated margins of error. Interpet the interval in context of the data. It may be helpful to create new data sets for each of the two countries first, and then use these data sets to construct the confidence intervals. Answer:

library(tidyverse)
library(infer)

set.seed(123)

# -------------------------------------------------------------------
# Variable 1: Physical Activity (Success = active 7 days a week)
# -------------------------------------------------------------------
active_data <- yrbss %>%
  filter(helmet_12m == "never") %>%
  drop_na(physically_active_7d) %>%
  mutate(active_ind = ifelse(physically_active_7d == 7, "yes", "no"))

ci_activity <- active_data %>%
  specify(response = active_ind, success = "yes") %>%
  generate(reps = 1000, type = "bootstrap") %>%
  calculate(stat = "prop") %>%
  get_ci(level = 0.95)

ci_activity

## # A tibble: 1 × 2
##   lower_ci upper_ci
##      <dbl>    <dbl>
## 1    0.293    0.314

moe_activity <- (ci_activity$upper_ci - ci_activity$lower_ci) / 2
moe_activity

## [1] 0.01054173

# -------------------------------------------------------------------
# Variable 2: Sleep Duration (Success = getting 8+ hours of sleep)
# Note: school_night_hours_sleep is a character variable with levels like "8", "9", "10+"
# -------------------------------------------------------------------
sleep_data <- yrbss %>%
  filter(helmet_12m == "never") %>%
  drop_na(school_night_hours_sleep) %>%
  mutate(sleep_ind = ifelse(school_night_hours_sleep %in% c("8", "9", "10+"), "yes", "no"))

ci_sleep <- sleep_data %>%
  specify(response = sleep_ind, success = "yes") %>%
  generate(reps = 1000, type = "bootstrap") %>%
  calculate(stat = "prop") %>%
  get_ci(level = 0.95)

ci_sleep

## # A tibble: 1 × 2
##   lower_ci upper_ci
##      <dbl>    <dbl>
## 1    0.291    0.312

moe_sleep <- (ci_sleep$upper_ci - ci_sleep$lower_ci) / 2
moe_sleep

## [1] 0.01074268

n <- 1000
p <- seq(from = 0, to = 1, by = 0.01)
me <- 2 * sqrt(p * (1 - p)/n)
dd <- data.frame(p = p, me = me)
ggplot(data = dd, aes(x = p, y = me)) + 
  geom_line() +
  labs(x = "Population Proportion", y = "Margin of Error")

Exercise 5

Describe the relationship between p and me. Include the margin of error vs. population proportion plot you constructed in your answer. For a given sample size, for which value of p is margin of error maximized? Answer: 1. Relationship Between \(p\) and Margin of Error (\(ME\)) The relationship between the population proportion \(p\) and the margin of error (\(ME\)) is symmetric, non-linear, and forms an inverted U-shaped parabolic curve.Key Observations:Symmetry around \(p = 0.5\): The curve is perfectly symmetrical about \(p = 0.5\). For instance, \(p = 0.2\) and \(p = 0.8\) produce the exact same margin of error because \(p(1 - p) = (0.2)(0.8) = (0.8)(0.2) = 0.16\).Minima at Extrema (\(p = 0\) and \(p = 1\)): At the extreme bounds where \(p = 0\) or \(p = 1\), the variability is zero (\(SE = 0\)), which makes the margin of error \(0\). If every individual in a population either completely possesses or lacks a trait, there is no sampling variability.Peak at \(p = 0.5\): As \(p\) moves toward \(0.5\) from either direction, the product \(p(1 - p)\) increases, reaching its absolute maximum at \(0.25\). Consequently, the margin of error reaches its maximum when \(p = 0.5\). 2. Margin of Error vs. Population Proportion Plot Below is the R code used to construct the dataset and visualization, followed by the resulting plot:

library(ggplot2)

# Set sample size
n <- 1000

# Generate sequence of proportions from 0 to 1
p <- seq(from = 0, to = 1, by = 0.01)

# Calculate margin of error
me <- 2 * sqrt(p * (1 - p) / n)

# Combine into data frame
dd <- data.frame(p = p, me = me)

# Create plot
ggplot(data = dd, aes(x = p, y = me)) +
  geom_line(color = "#2b5c8f", size = 1.2) +
  geom_point(aes(x = 0.5, y = 2 * sqrt(0.5 * 0.5 / 1000)), color = "red", size = 3) +
  annotate("text", x = 0.5, y = 0.033, label = "Maximum ME at p = 0.5 (ME = 0.0316)", color = "red", fontface = "bold") +
  labs(
    title = "Margin of Error vs. Population Proportion (n = 1000)",
    x = "Population Proportion (p)",
    y = "Margin of Error (ME)"
  ) +
  theme_minimal()

## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once per session.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

## Warning in geom_point(aes(x = 0.5, y = 2 * sqrt(0.5 * 0.5/1000)), color = "red", : All aesthetics have length 1, but the data has 101 rows.
## ℹ Please consider using `annotate()` or provide this layer with data containing
##   a single row.

3. Value of \(p\) Where Margin of Error is MaximizedFor a fixed sample size \(n\), the margin of error is maximized at \(p = 0.50\) (or \(50\%\)).Mathematical Explanation:The formula for the margin of error is:\[\text{ME} = 2 \times \sqrt{\frac{p(1 - p)}{n}}\]To find the value of \(p\) that maximizes \(\text{ME}\), we maximize the numerator function \(f(p) = p(1 - p) = p - p^2\). Taking the first derivative with respect to \(p\) and setting it to zero:\[\frac{df}{dp} = 1 - 2p = 0 \implies 2p = 1 \implies p = 0.50\]At \(n = 1000\) and \(p = 0.50\), the maximum margin of error is:\[\text{ME}_{\text{max}} = 2 \times \sqrt{\frac{0.5(1 - 0.5)}{1000}} = 2 \times \sqrt{\frac{0.25}{1000}} \approx 0.0316 \quad (3.16\%)\]

Exercise 6

Describe the sampling distribution of sample proportions at n=300 and p=0.1. Be sure to note the center, spread, and shape.

Answer:

To describe the sampling distribution of sample proportions (\(\hat{p}\)) for a sample size of \(n = 300\) and a population proportion of \(p = 0.10\), we evaluate the shape, center, and spread using the Central Limit Theorem for proportions.

1. Shape To check whether the sampling distribution of \(\hat{p}\) is approximately normal, we verify the Success-Failure Condition:\[n \cdot p = 300 \times 0.10 = 30 \ge 10\]\[n \cdot (1 - p) = 300 \times (1 - 0.10) = 300 \times 0.90 = 270 \ge 10\]Since both expected successes (\(30\)) and expected failures (\(270\)) are at least \(10\), the condition is satisfied.Shape: The sampling distribution of \(\hat{p}\) is unimodal and approximately symmetric (Normal): \(\hat{p} \sim N(\mu_{\hat{p}}, \sigma_{\hat{p}})\).

2. Center The mean (expected value) of the sampling distribution of sample proportions is equal to the true population proportion \(p\):\[\mu_{\hat{p}} = p = 0.10\]Center: The mean of the sampling distribution is centered at \(0.10\) (or \(10\%\)).

3. Spread (Standard Error)Assuming the sample of \(n = 300\) observations is selected randomly and represents less than \(10\%\) of the population (10% condition), the standard error (\(SE\)) measures the variability between sample proportions:\[\sigma_{\hat{p}} = SE = \sqrt{\frac{p(1 - p)}{n}}\]Calculation:\[\sigma_{\hat{p}} = \sqrt{\frac{0.10 \times (1 - 0.10)}{300}} = \sqrt{\frac{0.10 \times 0.90}{300}} = \sqrt{\frac{0.09}{300}} = \sqrt{0.0003} \approx 0.01732\]Spread: The standard error of the sampling distribution is approximately \(0.01732\) (or \(1.732\%\)).

# Insert code for Exercise 6 here

Exercise 7

Keep n constant and change p . How does the shape, center, and spread of the sampling distribution vary as p changes. You might want to adjust min and max for the x -axis for a better view of the distribution.

Answer: When holding the sample size (\(n\)) constant and varying the population proportion (\(p\)), the sampling distribution of sample proportions (\(\hat{p}\)) changes in shape, center, and spread in the following systematic ways: 1. Center (\(\mu_{\hat{p}}\))Behavior: The center of the sampling distribution always equals the population proportion:\[\mu_{\hat{p}} = p\]Effect of changing \(p\): As \(p\) increases, the peak of the distribution shifts directly to the right along the horizontal axis (\(x\)-axis). As \(p\) decreases, it shifts to the left.

2. Spread / Standard Error (\(\sigma_{\hat{p}}\))Behavior:

The variability (standard error) is calculated as:\[\sigma_{\hat{p}} = \sqrt{\frac{p(1 - p)}{n}}\]Effect of changing \(p\):Maximized at \(p = 0.50\): The standard error reaches its maximum value when \(p = 0.50\), where \(p(1 - p) = 0.25\). At \(p = 0.50\), the distribution has its widest spread.Decreases toward Extrema (\(p \to 0\) or \(p \to 1\)): As \(p\) moves further away from \(0.50\) toward either \(0\) or \(1\), \(p(1 - p)\) decreases, causing the standard error to shrink and making the distribution narrower and taller.Symmetry in Spread: Values of \(p\) equidistant from \(0.50\) (such as \(p = 0.20\) and \(p = 0.80\)) produce the exact same standard error and spread width.

3. Shape

Behavior: The shape is governed by the Success-Failure Condition (\(np \ge 10\) and \(n(1 - p) \ge 10\)).Effect of changing \(p\):Symmetric at \(p = 0.50\): When \(p = 0.50\), the distribution is perfectly symmetric and bell-shaped, regardless of \(n\) (so long as \(n \ge 20\)).Skewness near Extrema: As \(p\) gets close to \(0\) or \(1\), the distribution becomes skewed:When \(p\) is small (close to \(0\)), \(np < 10\) makes the distribution right-skewed (bounded by \(0\)).When \(p\) is large (close to \(1\)), \(n(1 - p) < 10\) makes the distribution left-skewed (bounded by \(1\)).Approximate Normality: As long as \(p\) stays within a moderate range such that both \(np \ge 10\) and \(n(1 - p) \ge 10\) hold, the distribution remains approximately normal.

library(ggplot2)

n <- 300
p_vals <- c(0.10, 0.50, 0.90)

# Generate simulated sampling distributions
set.seed(123)
sim_data <- map_dfr(p_vals, function(p) {
  tibble(
    p_label = paste0("p = ", p),
    phat = rbinom(10000, size = n, prob = p) / n
  )
})

# Plot showing shift in center, change in spread, and symmetry
ggplot(sim_data, aes(x = phat, fill = p_label)) +
  geom_density(alpha = 0.5) +
  scale_x_continuous(limits = c(0, 1), breaks = seq(0, 1, 0.1)) +
  labs(
    title = "Sampling Distributions of p-hat for n = 300 across varying p",
    x = "Sample Proportion (p-hat)",
    y = "Density",
    fill = "Parameter p"
  ) +
  theme_minimal()

Exercise 8

Now also change n . How does n appear to affect the distribution of p^ ? Answer: When keeping \(p\) constant and changing the sample size \(n\), \(n\) affects the sampling distribution of \(\hat{p}\) across shape, center, and spread in the following ways:

1. Spread (Standard Error \(\sigma_{\hat{p}}\)) — Inversely Proportional to $

\(Behavior: The standard error formula is:\)\(\sigma_{\hat{p}} = \sqrt{\frac{p(1 - p)}{n}}\)$Effect of Increasing \(n\):As \(n\) increases, the denominator gets larger, which decreases the standard error (\(\sigma_{\hat{p}}\)).The sampling distribution becomes narrower, taller, and more tightly clustered around the true population mean \(p\).To cut the standard error (and the margin of error) in half, you must increase the sample size \(n\) by a factor of 4.Effect of Decreasing \(n\):As \(n\) decreases, standard error increases, causing the distribution to spread out and become wider/flatter.

2. Shape — Approaches Normality as \(n\) IncreasesBehavior: The Central Limit Theorem dictates that larger sample sizes force the sampling distribution of \(\hat{p}\) to become more bell-shaped and symmetric.Success-Failure Condition (\(np \ge 10\) and \(n(1 - p) \ge 10\)):Small \(n\): If \(n\) is small, the expected counts \(np\) or \(n(1 - p)\) may fall below \(10\) (especially when \(p\) is close to \(0\) or \(1\)). This causes the distribution to be skewed (right-skewed if \(p < 0.5\), left-skewed if \(p > 0.5\)).Large \(n\): As \(n\) grows large enough to satisfy \(np \ge 10\) and \(n(1 - p) \ge 10\), any initial skewness disappears, and the shape becomes approximately Normal.3. Center (\(\mu_{\hat{p}}\)) — UnchangedBehavior: The mean (expected value) of the sampling distribution remains equal to the true population proportion regardless of sample size:\[\mu_{\hat{p}} = p\]Effect of changing \(n\): \(n\) has no effect on the center. The distribution will remain centered at \(p\) whether \(n = 30\) or \(n = 3,000\).

library(ggplot2)
library(dplyr)
library(purrr)

set.seed(123)
p <- 0.10
n_vals <- c(30, 100, 300, 1200)

# Generate sampling distributions for different n
sim_n_data <- map_dfr(n_vals, function(n_size) {
  tibble(
    n_label = factor(paste0("n = ", n_size), levels = paste0("n = ", n_vals)),
    phat = rbinom(10000, size = n_size, prob = p) / n_size
  )
})

# Overlay density curves
ggplot(sim_n_data, aes(x = phat, fill = n_label)) +
  geom_density(alpha = 0.4) +
  scale_x_continuous(limits = c(0, 0.3)) +
  labs(
    title = "Effect of Sample Size (n) on Sampling Distribution of p-hat (p = 0.10)",
    x = "Sample Proportion (p-hat)",
    y = "Density",
    fill = "Sample Size (n)"
  ) +
  theme_minimal()

## Warning: Removed 6 rows containing non-finite outside the scale range
## (`stat_density()`).

Exercise 9

Is there convincing evidence that those who sleep 10+ hours per day are more likely to strength train every day of the week? As always, write out the hypotheses for any tests you conduct and outline the status of the conditions for inference. If you find a significant difference, also quantify this difference with a confidence interval. Answer:

set.seed(123)

# 1. Prepare Data
yrbss_strength <- yrbss %>%
  drop_na(school_night_hours_sleep, strength_training_7d) %>%
  mutate(
    sleep_10plus = ifelse(school_night_hours_sleep == "10+", "10+ hours", "<10 hours"),
    sleep_10plus = factor(sleep_10plus, levels = c("10+ hours", "<10 hours")),
    strength_daily = ifelse(strength_training_7d == 7, "yes", "no"),
    strength_daily = factor(strength_daily, levels = c("yes", "no"))
  )

# 2. Compute Observed Difference in Proportions
obs_diff <- yrbss_strength %>%
  specify(strength_daily ~ sleep_10plus, success = "yes") %>%
  calculate(stat = "diff in props", order = c("10+ hours", "<10 hours"))

obs_diff

## Response: strength_daily (factor)
## Explanatory: sleep_10plus (factor)
## # A tibble: 1 × 1
##    stat
##   <dbl>
## 1 0.105

# 3. Conduct Hypothesis Test via Permutation / Simulation
null_dist <- yrbss_strength %>%
  specify(strength_daily ~ sleep_10plus, success = "yes") %>%
  hypothesize(null = "independence") %>%
  generate(reps = 1000, type = "permute") %>%
  calculate(stat = "diff in props", order = c("10+ hours", "<10 hours"))

# Calculate p-value
p_value <- null_dist %>%
  get_p_value(obs_stat = obs_diff, direction = "greater")

## Warning: Please be cautious in reporting a p-value of 0. This result is an approximation
## based on the number of `reps` chosen in the `generate()` step.
## ℹ See `get_p_value()` (`?infer::get_p_value()`) for more information.

p_value

## # A tibble: 1 × 1
##   p_value
##     <dbl>
## 1       0

# 4. Quantify Difference with a 95% Confidence Interval (Bootstrap)
boot_dist <- yrbss_strength %>%
  specify(strength_daily ~ sleep_10plus, success = "yes") %>%
  generate(reps = 1000, type = "bootstrap") %>%
  calculate(stat = "diff in props", order = c("10+ hours", "<10 hours"))

ci_diff <- boot_dist %>%
  get_ci(level = 0.95)

ci_diff

## # A tibble: 1 × 2
##   lower_ci upper_ci
##      <dbl>    <dbl>
## 1   0.0560    0.155

Exercise 10

Let’s say there has been no difference in likeliness to strength train every day of the week for those who sleep 10+ hours. What is the probablity that you could detect a change (at a significance level of 0.05) simply by chance? Hint: Review the definition of the Type 1 error.

Answer:

The probability of detecting a change simply by chance is equal to the significance level (\(\alpha\)), which is \(0.05\) (or \(5\%\)).

Explanation & Context:

Definition of Type I Error (\(\alpha\)):\[\text{Type I Error} = P(\text{Reject } H_0 \mid H_0 \text{ is true}) = \alpha\]The significance level \(\alpha\) represents the threshold set by the researcher for how much risk of a false positive they are willing to accept.

In Context:

When conducting a hypothesis test with \(\alpha = 0.05\), you are establishing a decision rule where \(5\%\) of all samples drawn from a population with no true difference will randomly produce a sample statistic extreme enough to fall in the rejection region, leading to a false conclusion of a difference.

Exercise 11

Suppose you’re hired by the local government to estimate the proportion of residents that attend a religious service on a weekly basis. According to the guidelines, the estimate must have a margin of error no greater than 1% with 95% confidence. You have no idea what to expect for p . How many people would you have to sample to ensure that you are within the guidelines? Hint: Refer to your plot of the relationship between p and margin of error. This question does not require using a dataset.

Answer:

To determine the required sample size (\(n\)) under the strictest conditions, we use the worst-case scenario for variability based on the relationship between population proportion (\(p\)) and margin of error (\(ME\)).

1. Determining the Value of \(p\) As established in the \(ME\) vs. \(p\) curve, the margin of error reaches its absolute maximum when \(p = 0.50\) (\(50\%\)).

When you have no prior information or expectation about \(p\), using \(p = 0.50\) provides the most conservative (worst-case) estimate. This guarantees that your sample size will be large enough to achieve a margin of error of \(1\%\) or less, regardless of what the true population proportion turns out to be.

2. Given Parameters

Target Margin of Error (\(ME\)): \(0.01\) (\(1\%\))

Confidence Level: \(95\%\)

Critical Value (\(z^*\)): \(1.96\) (for a \(95\%\) confidence level)

Assumed Proportion (\(p\)): \(0.50\)

3. Sample Size Formula & Calculation

The formula for the margin of error of a proportion is:\[ME = z^* \times \sqrt{\frac{p(1 - p)}{n}}\]Solving for sample size \(n\):\[n = \left( \frac{z^*}{ME} \right)^2 \times p(1 - p)\]Step-by-Step Substitution:\[n = \left( \frac{1.96}{0.01} \right)^2 \times (0.50)(1 - 0.50)\]\[n = (196)^2 \times 0.25\]\[n = 38,416 \times 0.25 = 9,604\]

Conclusion

It would need to sample at least 9,604 residents to ensure that your estimate of the proportion of residents attending weekly religious services has a margin of error no greater than \(1\%\) at a \(95\%\) confidence level.