set.seed(1234)
In this lab, you will investigate the ways in which the statistics from a random sample of data can serve as point estimates for population parameters. We’re interested in formulating a sampling distribution of our estimate in order to learn about the properties of the estimate, such as its distribution.
Setting a seed: We will take some random samples and build sampling distributions in this lab, which means you should set a seed at the start of your lab. If this concept is new to you, review the lab on probability.
In this lab, we will explore and visualize the data using the tidyverse suite of packages. We will also use the infer package for resampling.
Let’s load the packages.
library(tidyverse)
library(openintro)
library(infer)
A 2019 Gallup report states the following:
The premise that scientific progress benefits people has been embodied in discoveries throughout the ages – from the development of vaccinations to the explosion of technology in the past few decades, resulting in billions of supercomputers now resting in the hands and pockets of people worldwide. Still, not everyone around the world feels science benefits them personally.
The Wellcome Global Monitor finds that 20% of people globally do not believe that the work scientists do benefits people like them. In this lab, you will assume this 20% is a true population proportion and learn about how sample proportions can vary from sample to sample by taking smaller samples from the population. We will first create our population assuming a population size of 100,000. This means 20,000 (20%) of the population think the work scientists do does not benefit them personally and the remaining 80,000 think it does.
global_monitor <- tibble(
scientist_work = c(rep("Benefits", 80000), rep("Doesn't benefit", 20000))
)
The name of the data frame is global_monitor
and the
name of the variable that contains responses to the question “Do you
believe that the work scientists do benefit people like you?” is
scientist_work
.
We can quickly visualize the distribution of these responses using a bar plot.
ggplot(global_monitor, aes(x = scientist_work)) +
geom_bar() +
labs(
x = "", y = "",
title = "Do you believe that the work scientists do benefit people like you?"
) +
coord_flip()
We can also obtain summary statistics to confirm we constructed the data frame correctly.
global_monitor %>%
count(scientist_work) %>%
mutate(p = n /sum(n))
## # A tibble: 2 × 3
## scientist_work n p
## <chr> <int> <dbl>
## 1 Benefits 80000 0.8
## 2 Doesn't benefit 20000 0.2
In this lab, you have access to the entire population, but this is rarely the case in real life. Gathering information on an entire population is often extremely costly or impossible. Because of this, we often take a sample of the population and use that to understand the properties of the population.
If you are interested in estimating the proportion of people who
don’t think the work scientists do benefits them, you can use the
sample_n
command to survey the population.
samp1 <- global_monitor %>%
sample_n(50)
This command collects a simple random sample of size 50 from the
global_monitor
dataset, and assigns the result to
samp1
. This is similar to randomly drawing names from a hat
that contains the names of all in the population. Working with these 50
names is considerably simpler than working with all 100,000 people in
the population.
sample_n
function takes
a random sample of observations (i.e. rows) from the dataset, you can
still refer to the variables in the dataset with the same names. Code
you presented earlier for visualizing and summarizing the population
data will still be useful for the sample, however be careful to not
label your proportion p
since you’re now calculating a
sample statistic, not a population parameters. You can customize the
label of the statistics to indicate that it comes from the sample.#ANSWER #1
ggplot(samp1, aes(x = scientist_work)) +
geom_bar() +
labs(
x = "", y = "",
title = "SAMPLE OF 50 PEOPLE Do you believe that the work scientists do benefit people like you?"
) +
coord_flip()
samp1 %>%
count(scientist_work) %>%
mutate(pHAT = n /sum(n))
## # A tibble: 2 × 3
## scientist_work n pHAT
## <chr> <int> <dbl>
## 1 Benefits 37 0.74
## 2 Doesn't benefit 13 0.26
If you’re interested in estimating the proportion of all people who do not believe that the work scientists do benefits them, but you do not have access to the population data, your best single guess is the sample mean.
samp1 %>%
count(scientist_work) %>%
mutate(p_hat = n /sum(n))
## # A tibble: 2 × 3
## scientist_work n p_hat
## <chr> <int> <dbl>
## 1 Benefits 37 0.74
## 2 Doesn't benefit 13 0.26
Depending on which 50 people you selected, your estimate could be a bit above or a bit below the true population proportion of 0.26. In general, though, the sample proportion turns out to be a pretty good estimate of the true population proportion, and you were able to get it by sampling less than 1% of the population.
Would you expect the sample proportion to match the sample proportion of another student’s sample? Why, or why not? If the answer is no, would you expect the proportions to be somewhat different or very different? Ask a student team to confirm your answer. #The answer depends on the size of the sample, as the sample gets larger one would expect the differences between the two samples to become closer and closer. However, assuming a sample size of 50, the proportions will almost certainly be different but will stay relatively close to one another as the sample size is not extremely small. #ANSWER QUESTION 2
Take a second sample, also of size 50, and call it
samp2
. How does the sample proportion of samp2
compare with that of samp1
? Suppose we took two more
samples, one of size 100 and one of size 1000. Which would you think
would provide a more accurate estimate of the population proportion?
#ANSWER #3, the samples of 50 ( Samp1 / Samp3) differ only slightly but
the general proportion is roughly equivalent. If we took samples of size
100 and 1000, they would theoretically get more accurate and closer to
the 80/20 proportion.
set.seed(23)
samp3 <- global_monitor %>%
sample_n(50)
samp3 %>%
count(scientist_work) %>%
mutate(p_hat1 = n /sum(n))
## # A tibble: 2 × 3
## scientist_work n p_hat1
## <chr> <int> <dbl>
## 1 Benefits 39 0.78
## 2 Doesn't benefit 11 0.22
ggplot(samp3, aes(x = scientist_work)) +
geom_bar() +
labs(
x = "", y = "",
title = "SAMPLE OF 50 PEOPLE Do you believe that the work scientists do benefit people like you?"
) +
coord_flip()
Not surprisingly, every time you take another random sample, you might
get a different sample proportion. It’s useful to get a sense of just
how much variability you should expect when estimating the population
mean this way. The distribution of sample proportions, called the
sampling distribution (of the proportion), can help you
understand this variability. In this lab, because you have access to the
population, you can build up the sampling distribution for the sample
proportion by repeating the above steps many times. Here, we use R to
take 15,000 different samples of size 50 from the population, calculate
the proportion of responses in each sample, filter for only the
Doesn’t benefit responses, and store each result in a vector
called
sample_props50
. Note that we specify that
replace = TRUE
since sampling distributions are constructed
by sampling with replacement.
sample_props50 <- global_monitor %>%
rep_sample_n(size = 50, reps = 15000, replace = TRUE) %>%
count(scientist_work) %>%
mutate(p_hat = n /sum(n)) %>%
filter(scientist_work == "Doesn't benefit")
And we can visualize the distribution of these proportions with a histogram.
ggplot(data = sample_props50, aes(x = p_hat)) +
geom_histogram(binwidth = 0.02) +
labs(
x = "p_hat (Doesn't benefit)",
title = "Sampling distribution of p_hat",
subtitle = "Sample size = 50, Number of samples = 15000"
)
Next, you will review how this set of code works.
sample_props50
? Describe
the sampling distribution, and be sure to specifically note its center.
Make sure to include a plot of the distribution in your answer. #ANSWER
4 : There are two elements : people who don’t believe the scientist’s
work does/doesn’t benefit them. Center seems to be centered around 20%
of people not believing that scientists work benefits them, with some
margin of error extending to the right/left up to 40% and as low as 5%
but the vast majority of observations center around 20% which is also
the population proportion.ggplot(data = sample_props50, aes(x = p_hat)) +
geom_histogram(binwidth = 0.02) +
labs(
x = "p_hat (Doesn't benefit)",
title = "Sampling distribution of p_hat",
subtitle = "Sample size = 50, Number of samples = 15000"
)
The idea behind the rep_sample_n
function is
repetition. Earlier, you took a single sample of size
n
(50) from the population of all people in the population.
With this new function, you can repeat this sampling procedure
rep
times in order to build a distribution of a series of
sample statistics, which is called the sampling
distribution.
Note that in practice one rarely gets to build true sampling distributions, because one rarely has access to data from the entire population.
Without the rep_sample_n
function, this would be
painful. We would have to manually run the following code 15,000
times
global_monitor %>%
sample_n(size = 50, replace = TRUE) %>%
count(scientist_work) %>%
mutate(p_hat = n /sum(n)) %>%
filter(scientist_work == "Doesn't benefit")
## # A tibble: 1 × 3
## scientist_work n p_hat
## <chr> <int> <dbl>
## 1 Doesn't benefit 3 0.06
as well as store the resulting sample proportions each time in a separate vector.
Note that for each of the 15,000 times we computed a proportion, we did so from a different sample!
rep_sample_n
function does, try
modifying the code to create a sampling distribution of 25
sample proportions from samples of size 10,
and put them in a data frame named sample_props_small
.
Print the output. How many observations are there in this object called
sample_props_small
? What does each observation represent?
#ANSWER 5 there are 250 observations, each observation is one persons
opinion on whether scientific research benefits them or not ## Sample
size and the sampling distributionsample_props_small <- global_monitor %>%
rep_sample_n(size = 10, reps = 25, replace = TRUE) %>%
count(scientist_work) %>%
mutate(p_hat = n /sum(n)) %>%
filter(scientist_work == "Doesn't benefit")
ggplot(data = sample_props50, aes(x = p_hat)) +
geom_histogram(binwidth = 0.02) +
labs(
x = "p_hat (Doesn't benefit)",
title = "Sampling distribution of p_hat",
subtitle = "Sample size = 10, Number of samples = 25"
)
Mechanics aside, let’s return to the reason we used the
rep_sample_n
function: to compute a sampling distribution,
specifically, the sampling distribution of the proportions from samples
of 50 people.
ggplot(data = sample_props50, aes(x = p_hat)) +
geom_histogram(binwidth = 0.02)
The sampling distribution that you computed tells you much about estimating the true proportion of people who think that the work scientists do doesn’t benefit them. Because the sample proportion is an unbiased estimator, the sampling distribution is centered at the true population proportion, and the spread of the distribution indicates how much variability is incurred by sampling only 50 people at a time from the population.
In the remainder of this section, you will work on getting a sense of the effect that sample size has on your sampling distribution.
#haing trouble using the app, as mentioned above each observation is someone’s opinion on whether scientific research benefits them. As the sample size increases, the standard error tightens/becomes smaller, the mean of the sampling distribution goes closer to 20%, and the shape looks more and more like the standard normal distribution
sample_props_test <- global_monitor %>%
rep_sample_n(size = 10, reps = 5000, replace = TRUE) %>%
count(scientist_work) %>%
mutate(p_hat = n /sum(n)) %>%
filter(scientist_work == "Doesn't benefit")
ggplot(data = sample_props_test, aes(x = p_hat)) +
geom_histogram(binwidth = 0.02) +
labs(
x = "p_hat (Doesn't benefit)",
title = "Sampling distribution of p_hat",
subtitle = "Sample size = 10, Number of samples = 5000"
)
#Just change the 10 to a 50 and 100 in order to increase the size of the samples
So far, you have only focused on estimating the proportion of those you think the work scientists doesn’t benefit them. Now, you’ll try to estimate the proportion of those who think it does.
Note that while you might be able to answer some of these questions using the app, you are expected to write the required code and produce the necessary plots and summary statistics. You are welcome to use the app for exploration.
sample_props_test_7 <- global_monitor %>%
rep_sample_n(size = 15, reps = 1, replace = TRUE) %>%
count(scientist_work) %>%
mutate(p_hat = n /sum(n)) %>%
filter(scientist_work == "Doesn't benefit")
ggplot(data = sample_props_test_7, aes(x = p_hat)) +
geom_histogram(binwidth = 0.02) +
labs(
x = "p_hat (Doesn't benefit)",
title = "Sampling distribution of p_hat",
subtitle = "Sample size = 15, Number of samples = 1"
)
sample_props15
. Plot the data, then
describe the shape of this sampling distribution. Based on this sampling
distribution, what would you guess the true proportion of those who
think the work scientists do enchances their lives to be? Finally,
calculate and report the population proportion.sample_props15 <- global_monitor %>%
rep_sample_n(size = 15, reps = 2000, replace = TRUE) %>%
count(scientist_work) %>%
mutate(p_hat = n /sum(n)) %>%
filter(scientist_work == "Benefits")
ggplot(data = sample_props15, aes(x = p_hat)) +
geom_histogram(binwidth = 0.01) +
labs(
x = "p_hat (Doesn't benefit)",
title = "Sampling distribution of p_hat",
subtitle = "Sample size = 15, Number of samples = 2000"
)
sample_props150
. Describe the shape
of this sampling distribution and compare it to the sampling
distribution for a sample size of 15. Based on this sampling
distribution, what would you guess to be the true proportion of those
who think the work scientists do enchances their lives?#Answer 9 : The sampling distribution is slighlty more symmetrical, although it still has some skew. Seems like the true proportion is 80% for people who believe science benefits. As the sample size of the repeated samples increase, seems like the distribution becomes more and more normal.
sample_props150 <- global_monitor %>%
rep_sample_n(size = 150, reps = 2000, replace = TRUE) %>%
count(scientist_work) %>%
mutate(p_hat = n /sum(n)) %>%
filter(scientist_work == "Benefits")
ggplot(data = sample_props15, aes(x = p_hat)) +
geom_histogram(binwidth = 0.01) +
labs(
x = "p_hat (Doesn't benefit)",
title = "Sampling distribution of p_hat",
subtitle = "Sample size = 150, Number of samples = 2000"
)