Getting Started

Load packages

In this lab, we will explore and visualize the data using the tidyverse suite of packages, and perform statistical inference using infer. The data can be found in the companion package for OpenIntro resources, openintro.

Let’s load the packages.

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

The data

Every two years, the Centers for Disease Control and Prevention conduct the Youth Risk Behavior Surveillance System (YRBSS) survey, where it takes data from high schoolers (9th through 12th grade), to analyze health patterns. You will work with a selected group of variables from a random sample of observations during one of the years the YRBSS was conducted.

Load the yrbss data set into your workspace.

data('yrbss', package='openintro')

There are observations on 13 different variables, some categorical and some numerical. The meaning of each variable can be found by bringing up the help file:

?yrbss
  1. What are the cases in this data set? How many cases are there in our sample?

There are 13 cases which describe different health attributes and habits and 13,583 rows for each case

Remember that you can answer this question by viewing the data in the data viewer or by using the following command:

glimpse(yrbss)
## Rows: 13,583
## Columns: 13
## $ age                      <int> 14, 14, 15, 15, 15, 15, 15, 14, 15, 15, 15, 1…
## $ gender                   <chr> "female", "female", "female", "female", "fema…
## $ grade                    <chr> "9", "9", "9", "9", "9", "9", "9", "9", "9", …
## $ hispanic                 <chr> "not", "not", "hispanic", "not", "not", "not"…
## $ race                     <chr> "Black or African American", "Black or Africa…
## $ height                   <dbl> NA, NA, 1.73, 1.60, 1.50, 1.57, 1.65, 1.88, 1…
## $ weight                   <dbl> NA, NA, 84.37, 55.79, 46.72, 67.13, 131.54, 7…
## $ helmet_12m               <chr> "never", "never", "never", "never", "did not …
## $ text_while_driving_30d   <chr> "0", NA, "30", "0", "did not drive", "did not…
## $ physically_active_7d     <int> 4, 2, 7, 0, 2, 1, 4, 4, 5, 0, 0, 0, 4, 7, 7, …
## $ hours_tv_per_school_day  <chr> "5+", "5+", "5+", "2", "3", "5+", "5+", "5+",…
## $ strength_training_7d     <int> 0, 0, 0, 0, 1, 0, 2, 0, 3, 0, 3, 0, 0, 7, 7, …
## $ school_night_hours_sleep <chr> "8", "6", "<5", "6", "9", "8", "9", "6", "<5"…

Exploratory data analysis

You will first start with analyzing the weight of the participants in kilograms: weight.

Using visualization and summary statistics, describe the distribution of weights. The summary function can be useful.

summary(yrbss$weight)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
##   29.94   56.25   64.41   67.91   76.20  180.99    1004
  1. How many observations are we missing weights from?

There are 1004 observations missing form the weights column

sum(is.na(yrbss$weight))
## [1] 1004

Next, consider the possible relationship between a high schooler’s weight and their physical activity. Plotting the data is a useful first step because it helps us quickly visualize trends, identify strong associations, and develop research questions.

First, let’s create a new variable physical_3plus, which will be coded as either “yes” if they are physically active for at least 3 days a week, and “no” if not.

yrbss <- yrbss %>% 
  mutate(physical_3plus = ifelse(yrbss$physically_active_7d > 2, "yes", "no"))
yrbss2 <- na.omit(yrbss)
  1. Make a side-by-side boxplot of physical_3plus and weight. Is there a relationship between these two variables? What did you expect and why?

I expected that those who were more active would weigh less than those who are not physically active because when you’re active, you burn calories, and therefore lose weight, but I guess there are other factors that’s not just about how active you are but also about how much you eat. The box plot shows that those who are less active weigh less.

yrbss2 <- yrbss %>% 
  mutate(physical_3plus = ifelse(yrbss$physically_active_7d > 2, "yes", "no")) %>%
  na.exclude()
ggplot(yrbss2, aes(x=weight, y=physical_3plus)) + geom_boxplot() + theme_bw()

The box plots show how the medians of the two distributions compare, but we can also compare the means of the distributions using the following to first group the data by the physical_3plus variable, and then calculate the mean weight in these groups using the mean function while ignoring missing values by setting the na.rm argument to TRUE.

yrbss2 %>%
  group_by(physical_3plus) %>%
  summarise(mean_weight = mean(weight, na.rm = TRUE))
## # A tibble: 2 × 2
##   physical_3plus mean_weight
##   <chr>                <dbl>
## 1 no                    67.1
## 2 yes                   68.7

There is an observed difference, but is this difference statistically significant? In order to answer this question we will conduct a hypothesis test.

Inference

  1. Are all conditions necessary for inference satisfied? Comment on each. You can compute the group sizes with the summarize command above by defining a new variable with the definition n().

Both independence and normality have been satisfied because the data is a sample of a random distribution of students from across different schools which shows independence, and based on the box plots distribution, we can assume that the data is normal

 yrbss2 %>%
  group_by(physical_3plus) %>%
  summarise(mean_weight = mean(weight, na.rm = TRUE), count = n())
## # A tibble: 2 × 3
##   physical_3plus mean_weight count
##   <chr>                <dbl> <int>
## 1 no                    67.1  2656
## 2 yes                   68.7  5695
  1. Write the hypotheses for testing if the average weights are different for those who exercise at least times a week and those who don’t.

H0: Students who engage in physical activity for at least 3 days a week do not exhibit a distinct difference in average weight when compared to their counterparts who are not as physically active for at least 3 days a week.. HA:Students who are physically active for at least 3 days a week have a different average weight compared to those who are not as physically active for the same duration.

Next, we will introduce a new function, hypothesize, that falls into the infer workflow. You will use this method for conducting hypothesis tests.

But first, we need to initialize the test, which we will save as obs_diff.

obs_diff <- yrbss2 %>%
  drop_na(physical_3plus) %>%
  specify(weight ~ physical_3plus) %>%
  calculate(stat = "diff in means", order = c("yes", "no"))

Notice how you can use the functions specify and calculate again like you did for calculating confidence intervals. Here, though, the statistic you are searching for is the difference in means, with the order being yes - no != 0.

After you have initialized the test, you need to simulate the test on the null distribution, which we will save as null.

null_dist <- yrbss2 %>%
  drop_na(physical_3plus) %>%
  specify(weight ~ physical_3plus) %>%
  hypothesize(null = "independence") %>%
  generate(reps = 1000, type = "permute") %>%
  calculate(stat = "diff in means", order = c("yes", "no"))

Here, hypothesize is used to set the null hypothesis as a test for independence. In one sample cases, the null argument can be set to “point” to test a hypothesis relative to a point estimate.

Also, note that the type argument within generate is set to permute, which is the argument when generating a null distribution for a hypothesis test.

We can visualize this null distribution with the following code:

ggplot(data = null_dist, aes(x = stat)) +
  geom_histogram()

  1. How many of these null permutations have a difference of at least obs_stat?

None of them, 0

total_num <- sum(null_dist$stat >= obs_diff$stat)
print(total_num)
## [1] 0

Now that the test is initialized and the null distribution formed, you can calculate the p-value for your hypothesis test using the function get_p_value.

null_dist %>%
  get_p_value(obs_stat = obs_diff, direction = "two_sided")
## # A tibble: 1 × 1
##   p_value
##     <dbl>
## 1       0

This the standard workflow for performing hypothesis tests.

  1. Construct and record a confidence interval for the difference between the weights of those who exercise at least three times a week and those who don’t, and interpret this interval in context of the data.

With a 95% level of confidence, we can state that students who engage in physical activity at least three times a week have an average weight within the range of 68.09 kg to 68.8 kg. Likewise, students who do not participate in physical activity at least three times a week exhibit an average weight in the range of 66.13 kg to 67.22 kg with the same level of confidence.

## get sd for both values, 17.6 for no and 16.5 for yes
yrbss %>% 
  group_by(physical_3plus) %>% 
  summarise(sd_weight = sd(weight, na.rm = TRUE))
## # A tibble: 3 × 2
##   physical_3plus sd_weight
##   <chr>              <dbl>
## 1 no                  17.6
## 2 yes                 16.5
## 3 <NA>                17.6
yrbss %>% 
  group_by(physical_3plus) %>% 
  summarise(mean_weight = mean(weight, na.rm = TRUE))
## # A tibble: 3 × 2
##   physical_3plus mean_weight
##   <chr>                <dbl>
## 1 no                    66.7
## 2 yes                   68.4
## 3 <NA>                  69.9
yrbss %>% 
  group_by(physical_3plus) %>% 
  summarise(freq = table(weight)) %>%
  summarise(n = sum(freq))
## # A tibble: 3 × 2
##   physical_3plus     n
##   <chr>          <int>
## 1 no              4022
## 2 yes             8342
## 3 <NA>             215
mean_no <- 66.67389
size_no <- 4022
sd_no <- 17.63805
mean_yes <- 68.44847
size_yes <- 8342
sd_yes <- 16.47832

z = 1.96

up_no <- mean_no + z*(sd_no/sqrt(size_no))
low_not <- mean_no - z*(sd_no/sqrt(size_no))
up_no
## [1] 67.219
low_not
## [1] 66.12878
u_ci <- mean_yes + z*(sd_yes/sqrt(size_yes))
l_ci <- mean_yes - z*(sd_yes/sqrt(size_yes))

u_ci
## [1] 68.80209
l_ci
## [1] 68.09485

More Practice

  1. Calculate a 95% confidence interval for the average height in meters (height) and interpret it in context.

With 95% confident that the average height of the students in this population is between 1.689m and 1.693m.

height_table <- as.data.frame(table(yrbss$height))
total_frequency <- sum(height_table$Freq)

sample_mean_height <- mean(yrbss$height, na.rm = TRUE)
sample_sd_height <- sd(yrbss$height, na.rm = TRUE)

n_height <- yrbss %>% 
  summarise(freq = table(height)) %>%
  summarise(n = sum(freq, na.rm = TRUE))

z_value <- 1.96  
upper_ci_height <- sample_mean_height + z_value * (sample_sd_height / sqrt(n_height))
lower_ci_height <- sample_mean_height - z_value * (sample_sd_height / sqrt(n_height))

upper_ci_height
##          n
## 1 1.693071
  1. Calculate a new confidence interval for the same parameter at the 90% confidence level. Comment on the width of this interval versus the one obtained in the previous exercise.

The new range is 1.689705 to 1.692777. When we compare it to the first range of 1.689411 to 1.693071, we notice a small difference between them.

crit_90 <- 1.645
up_conf <- sample_mean_height + crit_90*(sample_sd_height/sqrt(n_height))
low_conf <- sample_mean_height - crit_90*(sample_sd_height/sqrt(n_height))
up_conf
##          n
## 1 1.692777
low_conf
##          n
## 1 1.689705
range_95 <- (upper_ci_height - lower_ci_height)
range_90 <- (up_conf - low_conf)
range_95
##             n
## 1 0.003659302
range_90
##           n
## 1 0.0030712
  1. Conduct a hypothesis test evaluating whether the average height is different for those who exercise at least three times a week and those who don’t.

The p-value of 0.05 is less than the significance level alpha of 0.1. Therefore, we reject the null hypothesis. This suggests a significant difference in the average height between individuals who engage in physical activity at least three times a week and those who do not. This finding may appear unexpected.

yrbss_no_active <- filter(yrbss2, physical_3plus == "no" & height != 0)
height_no_active <- yrbss_no_active$height

yrbss_active <- filter(yrbss2, physical_3plus == "yes" & weight != 0)
height_active <- yrbss_active$height

sample_size_active <- 5695
sample_size_no_active <- 2656

degrees_of_freedom <- sample_size_active - 1

mean_height_no_active <- mean(height_no_active)
mean_height_active <- mean(height_active)
sample_sd_heighteight_no_active <- sd(height_no_active)
sample_sd_heighteight_active <- sd(height_active)

standard_error <- sqrt((sample_sd_heighteight_active^2) / sample_size_active + (sample_sd_heighteight_no_active^2) / sample_size_no_active)

t_value <- qt(0.05/2, degrees_of_freedom, lower.tail = FALSE)

point_estimate <- mean_height_active - mean_height_no_active

lower_confidence_limit <- point_estimate - t_value * standard_error
upper_confidence_limit <- point_estimate + t_value * standard_error

lower_confidence_limit
## [1] 0.03329348
  1. Now, a non-inference task: Determine the number of different options there are in the dataset for the hours_tv_per_school_day there are.

There are 7 different options excluding the NA

yrbss %>%group_by(hours_tv_per_school_day)%>% summarise(n())
## # A tibble: 8 × 2
##   hours_tv_per_school_day `n()`
##   <chr>                   <int>
## 1 1                        1750
## 2 2                        2705
## 3 3                        2139
## 4 4                        1048
## 5 5+                       1595
## 6 <1                       2168
## 7 do not watch             1840
## 8 <NA>                      338
  1. Come up with a research question evaluating the relationship between height or weight and sleep. Formulate the question in a way that it can be answered using a hypothesis test and/or a confidence interval. Report the statistical results, and also provide an explanation in plain language. Be sure to check all assumptions, state your \(\alpha\) level, and conclude in context.

** Is there a significant difference in the average height between students who report getting enough sleep (more than 8 hours per night) and those who do not?

H0 (Null Hypothesis): The average height is the same for students who get enough sleep and those who do not.

HA (Alternative Hypothesis): The average height is different for students who get enough sleep compared to those who do not.

To test this hypothesis, we will conduct a two-sample t-test to compare the means of height for the two groups. We will set our significance level (\(\alpha\)) at 0.05.

Assumptions:

Independence: We assume that the data for height and sleep is independent, meaning one student’s height or sleep duration does not affect another student’s height or sleep duration. Normality: We assume that the heights in both groups (students who get enough sleep and those who do not) are normally distributed. **

enough_sleep <- yrbss %>% filter(school_night_hours_sleep > 8)  
not_enough_sleep <- yrbss %>% filter(school_night_hours_sleep <= 8)

t_test_result <- t.test(enough_sleep$height, not_enough_sleep$height, var.equal = TRUE)

t_test_result
## 
##  Two Sample t-test
## 
## data:  enough_sleep$height and not_enough_sleep$height
## t = -0.78715, df = 11479, p-value = 0.4312
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  -0.011173797  0.004770861
## sample estimates:
## mean of x mean of y 
##  1.687716  1.690918

Conclusion

since the p-value is less that the alpha, we reject the null hypothesis. This means that there is a significant difference in the average height between students who get enough sleep and people who do not