library(tidyverse)
library(openintro)
library(infer)
library(kableExtra)R-Lab - Inference for Categorical Data
Inference for Categorical Data
Getting Started
Load packages
Data
data("yrbss")
head(yrbss)# A tibble: 6 × 13
age gender grade hispanic race height weight helmet_12m
<int> <chr> <chr> <chr> <chr> <dbl> <dbl> <chr>
1 14 female 9 not Black or African American NA NA never
2 14 female 9 not Black or African American NA NA never
3 15 female 9 hispanic Native Hawaiian or Other… 1.73 84.4 never
4 15 female 9 not Black or African American 1.6 55.8 never
5 15 female 9 not Black or African American 1.5 46.7 did not r…
6 15 female 9 not Black or African American 1.57 67.1 did not r…
# ℹ 5 more variables: text_while_driving_30d <chr>, physically_active_7d <int>,
# hours_tv_per_school_day <chr>, strength_training_7d <int>,
# school_night_hours_sleep <chr>
Set Seed
set.seed(64488)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: Please see the lovely table below.
counts <- yrbss |>
count(text_while_driving_30d)
kbl(counts) |>
kable_styling() |>
row_spec(row = 0, color = "dodgerblue")| text_while_driving_30d | n |
|---|---|
| 0 | 4792 |
| 1-2 | 925 |
| 10-19 | 373 |
| 20-29 | 298 |
| 3-5 | 493 |
| 30 | 827 |
| 6-9 | 311 |
| did not drive | 4646 |
| NA | 918 |
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: The proportion of people who have texted while driving every day in the past 30 days and never wear helmets is 0.071.
no_helmet <- yrbss|>
filter(helmet_12m == "never") |>
mutate(text_ind =
ifelse(text_while_driving_30d == "30","yes","no"))
head(no_helmet)# A tibble: 6 × 14
age gender grade hispanic race height weight helmet_12m
<int> <chr> <chr> <chr> <chr> <dbl> <dbl> <chr>
1 14 female 9 not Black or African American NA NA never
2 14 female 9 not Black or African American NA NA never
3 15 female 9 hispanic Native Hawaiian or Other… 1.73 84.4 never
4 15 female 9 not Black or African American 1.6 55.8 never
5 14 male 9 not Black or African American 1.88 71.2 never
6 15 male 9 not Black or African American 1.75 63.5 never
# ℹ 6 more variables: text_while_driving_30d <chr>, physically_active_7d <int>,
# hours_tv_per_school_day <chr>, strength_training_7d <int>,
# school_night_hours_sleep <chr>, text_ind <chr>
no_helmet |>
filter(!is.na(text_ind)) |>
summarise(prop_text_no_helmet = mean(text_ind == "yes")) |>
pull()[1] 0.07119791
no_helmet <- no_helmet[!is.na(no_helmet$text_ind), ]Inference on proportions
set.seed(64488)
boot_n_h <- no_helmet |>
specify(response = text_ind, success = "yes") |>
generate(reps = 1000, type = "bootstrap") |>
calculate(stat = "prop")
get_ci(boot_n_h)# A tibble: 1 × 2
lower_ci upper_ci
<dbl> <dbl>
1 0.0649 0.0772
# margin_of_error and margin_of_error2 are really close in values and i understand the method in this chunk
ci <- get_ci(boot_n_h, level = 0.95)
margin_of_error <- (ci$upper_ci - ci$lower_ci) / 2
print(margin_of_error)[1] 0.006152929
# coding while unsure
standard_error2 <- sd(boot_n_h$stat)
margin_of_error2 <- 1.96 * standard_error2
margin_of_error2[1] 0.006065209
se <- boot_n_h |>
summarize(se = sd(stat)) |>
pull()
se[1] 0.003094494
me_se <- 1.96*se
me_se[1] 0.006065209
Exercise 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 for the proportion estimate is plus or minus 0.006155 (a 95% confidence interval of [0.06489, 0.07720])
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: In the chunks below, I checked for the confidence interval of those identified as males who texted and those who identified as hispanic and texted. As suggested in the lab directions, I created new data frames for each scenario before running calculations.
The margin of error for the male/text proportion is plus or minus 0.0085 (a 95% confidence interval of [0.0689, 0.0859])
The margin of error for the hispanic/text proportion is plus or minus 0.0112 (a 95% confidence interval of [0.0493, 0.0717])
xy_text <- no_helmet |>
filter(gender == "male")
hispanic_text <- no_helmet |>
filter(hispanic == "hispanic")xy_text |>
summarise(prop_xy_text = mean(text_ind == "yes")) |>
pull()[1] 0.07710652
set.seed(64488)
boot_xy_text <- xy_text |>
specify(response = text_ind, success = "yes") |>
generate(reps = 1000, type = "bootstrap") |>
calculate(stat = "prop")
ci <- get_ci(boot_xy_text, level = 0.95)
margin_of_error3 <- (ci$upper_ci - ci$lower_ci) / 2
print(margin_of_error3)[1] 0.008485692
hispanic_text |>
summarise(prop_h_text = mean(text_ind == "yes")) |>
pull()[1] 0.06077982
set.seed(64488)
boot_h_t <- hispanic_text |>
specify(response = text_ind, success = "yes") |>
generate(reps = 1000, type = "bootstrap") |>
calculate(stat = "prop")
ci <- get_ci(boot_h_t, level = 0.95)
margin_of_error4 <- (ci$upper_ci - ci$lower_ci) / 2
print(margin_of_error4)[1] 0.01118836
How does the proportion affect the margin of error?
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: For a given sample size, the margin of error is maximized when p and 1-p are equal. As the two values (p & 1-p) diverge, the margin of error becomes smaller. See above plot.
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 𝑝. 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 𝑝 and margin of error. This question does not require using a dataset.
Answer: Assuming p-hat is 0.50, we plug and chug using the formula in the chunk below. The answer is 9,604.
Out of curiosity, I wondered what would happen if we doubled the ME. The necessary sample size dropped by 75%!
n <- (0.5*0.5*1.96^2)/0.01^2
n[1] 9604
# what happens if we double the margin of error but keep everything else the same?
n_2 <- (0.5*0.5*1.96^2)/0.02^2
n_2[1] 2401
# what happens if we make p-hat closer to 0.05 based on the plot in Exercise 5? It seemed that when the proportion was around 0.05/0.95 it intersected the ME at 0.01.
n3 <- (0.05*0.95*1.96^2)/0.01^2
n3[1] 1824.76
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: 316 students sleep 10+ hours daily. The strength training variable is all integers.
The null hypothesis is that the proportion of (p-daily-gym) heavy sleepers who weight lift daily is the same as the proportion of (p-some-gym) heavy sleepers who don’t weight lift daily.
The alternative hypothesis is that p-daily-gym is greater than p-some-gym.
[gnashing of the teeth] Wouldn’t it make more sense to compare gym-7-days to gym-0, gym-1, gym-2-day… and so forth? Or to compare gym-7 to gym-0 to gym-some? I leave these notes as a testament to my effort.
heavy_sleep <- yrbss |>
filter(school_night_hours_sleep == "10+") |>
mutate(gym_time = ifelse(strength_training_7d == 7, "yes", "no")) |>
filter(!is.na(gym_time))
head(heavy_sleep)# A tibble: 6 × 14
age gender grade hispanic race height weight helmet_12m
<int> <chr> <chr> <chr> <chr> <dbl> <dbl> <chr>
1 15 female 9 not Black or African American 1.68 69.8 did not r…
2 18 male 12 not Black or African American 1.8 107. always
3 18 female 11 not Black or African American NA NA did not r…
4 18 male 11 not Black or African American 1.63 54.9 never
5 17 female 11 not Black or African American 1.65 62.6 did not r…
6 17 female 10 not Black or African American 1.65 61.2 never
# ℹ 6 more variables: text_while_driving_30d <chr>, physically_active_7d <int>,
# hours_tv_per_school_day <chr>, strength_training_7d <int>,
# school_night_hours_sleep <chr>, gym_time <chr>
set.seed(64488)
boot_g_y <- heavy_sleep |>
specify(response = gym_time, success = "yes") |>
generate(reps = 1000, type = "bootstrap") |>
calculate(stat = "prop")
get_ci(boot_g_y)# A tibble: 1 × 2
lower_ci upper_ci
<dbl> <dbl>
1 0.221 0.317
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: I wish to protest this question.