Foundations for statistical inference - Sampling distributions

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.

The data

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.

Source: World Science Day: Is Knowledge Power?

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

The unknown sampling distribution

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)

Exercise 1

Describe the distribution of responses in this sample. How does it compare to the distribution of responses in the population. Hint: Although the 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.

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           42  0.84
## 2 Doesn't benefit     8  0.16

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.

Exercise 2

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.

I would not expect the sample proportion to match the sample proportion of another student’s sample because sample proportion is denoted as a random variable.Each time the samplig is done the, outcome will almost always be different.

Exercise 3

Take a second sample, also of size 50, and call it samp2.

How does the sample proportion of samp2 compare with that of samp1?

Answer

The sample proportion of samp2 is similiar to samp1 because both sample sizes are the same and the ratios are close.

samp2 <- global_monitor %>%
sample_n(50)

samp2 %>%
  count(scientist_work) %>%
  mutate(p_hat = n /sum(n))
## # A tibble: 2 × 3
##   scientist_work      n p_hat
##   <chr>           <int> <dbl>
## 1 Benefits           44  0.88
## 2 Doesn't benefit     6  0.12

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

The sample with size 1000 provided a more accurate estimate of the population propotion because larger sample size provides a more accurate n and p_hat values.

samp3 <- global_monitor %>%
sample_n(100)

samp3 %>%
  count(scientist_work) %>%
  mutate(p_hat = n /sum(n))
## # A tibble: 2 × 3
##   scientist_work      n p_hat
##   <chr>           <int> <dbl>
## 1 Benefits           79  0.79
## 2 Doesn't benefit    21  0.21
samp4 <- global_monitor %>%
sample_n(1000)

samp4 %>%
  count(scientist_work) %>%
  mutate(p_hat = n /sum(n))
## # A tibble: 2 × 3
##   scientist_work      n p_hat
##   <chr>           <int> <dbl>
## 1 Benefits          781 0.781
## 2 Doesn't benefit   219 0.219

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.

Exercise 4

How many elements are there in sample_props50?

Answer

There are 15000 elements in sample_props50 because.

count(sample_props50)
## # A tibble: 15,000 × 2
## # Groups:   replicate [15,000]
##    replicate     n
##        <int> <int>
##  1         1     1
##  2         2     1
##  3         3     1
##  4         4     1
##  5         5     1
##  6         6     1
##  7         7     1
##  8         8     1
##  9         9     1
## 10        10     1
## # … with 14,990 more rows

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

The sampling distribution shows a distribution that is symmetric and with a center about 20%.

ggplot(data = sample_props50, aes(x = p_hat)) +
  geom_histogram(binwidth = 0.02) +
  labs(
    x = "p_hat (Doesn't benefit)",
    title = "Sampling distribution of sample_props50",
    subtitle = "Sample size = 50, Number of samples = 15000"
  )

Interlude: Sampling distributions

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    17  0.34

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!

Exercise 5

To make sure you understand how sampling distributions are built, and exactly what the 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.

Answer

How many observations are there in this object called sample_props_small? There are 25 observations in sample_props_small.

What does each observation represent? Each observation represents the sample proportion that makes up the sampling distribution(responses in each sample for population doesn’t believe in scientists).

sample_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")

sample_props_small
## # A tibble: 21 × 4
## # Groups:   replicate [21]
##    replicate scientist_work      n p_hat
##        <int> <chr>           <int> <dbl>
##  1         1 Doesn't benefit     2   0.2
##  2         2 Doesn't benefit     3   0.3
##  3         3 Doesn't benefit     2   0.2
##  4         5 Doesn't benefit     2   0.2
##  5         6 Doesn't benefit     2   0.2
##  6         7 Doesn't benefit     5   0.5
##  7         8 Doesn't benefit     2   0.2
##  8        10 Doesn't benefit     3   0.3
##  9        11 Doesn't benefit     2   0.2
## 10        12 Doesn't benefit     5   0.5
## # … with 11 more rows

Sample size and the sampling distribution

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)

Exercise 6

Use the app below to create sampling distributions of proportions of Doesn’t benefit from samples of size 10, 50, and 100. Use 5,000 simulations.

sampl_p10 <- 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")
sampl_p10
## # A tibble: 4,497 × 4
## # Groups:   replicate [4,497]
##    replicate scientist_work      n p_hat
##        <int> <chr>           <int> <dbl>
##  1         1 Doesn't benefit     2   0.2
##  2         2 Doesn't benefit     3   0.3
##  3         3 Doesn't benefit     2   0.2
##  4         4 Doesn't benefit     3   0.3
##  5         5 Doesn't benefit     2   0.2
##  6         6 Doesn't benefit     2   0.2
##  7         7 Doesn't benefit     3   0.3
##  8         9 Doesn't benefit     1   0.1
##  9        10 Doesn't benefit     1   0.1
## 10        11 Doesn't benefit     2   0.2
## # … with 4,487 more rows
mean(sampl_p10$p_hat)
## [1] 0.2238159
sd(sampl_p10$p_hat)
## [1] 0.1119618
sampl_p50 <- global_monitor %>%
                    rep_sample_n(size = 50, reps = 5000, replace = TRUE) %>%
   count(scientist_work) %>%
   mutate(p_hat = n /sum(n)) %>%
   filter(scientist_work == "Doesn't benefit")
sampl_p50
## # A tibble: 5,000 × 4
## # Groups:   replicate [5,000]
##    replicate scientist_work      n p_hat
##        <int> <chr>           <int> <dbl>
##  1         1 Doesn't benefit     8  0.16
##  2         2 Doesn't benefit     8  0.16
##  3         3 Doesn't benefit     8  0.16
##  4         4 Doesn't benefit    11  0.22
##  5         5 Doesn't benefit    14  0.28
##  6         6 Doesn't benefit     6  0.12
##  7         7 Doesn't benefit    10  0.2 
##  8         8 Doesn't benefit     9  0.18
##  9         9 Doesn't benefit    14  0.28
## 10        10 Doesn't benefit     7  0.14
## # … with 4,990 more rows
mean(sampl_p50$p_hat)
## [1] 0.199292
sd(sampl_p50$p_hat)
## [1] 0.05569667
sampl_p100 <- global_monitor %>%
                    rep_sample_n(size = 100, reps = 5000, replace = TRUE) %>%
   count(scientist_work) %>%
   mutate(p_hat = n /sum(n)) %>%
   filter(scientist_work == "Doesn't benefit")
sampl_p100
## # A tibble: 5,000 × 4
## # Groups:   replicate [5,000]
##    replicate scientist_work      n p_hat
##        <int> <chr>           <int> <dbl>
##  1         1 Doesn't benefit    18  0.18
##  2         2 Doesn't benefit    20  0.2 
##  3         3 Doesn't benefit    20  0.2 
##  4         4 Doesn't benefit    28  0.28
##  5         5 Doesn't benefit    18  0.18
##  6         6 Doesn't benefit    16  0.16
##  7         7 Doesn't benefit    25  0.25
##  8         8 Doesn't benefit    20  0.2 
##  9         9 Doesn't benefit    23  0.23
## 10        10 Doesn't benefit    17  0.17
## # … with 4,990 more rows
mean(sampl_p100$p_hat)
## [1] 0.199436
sd(sampl_p100$p_hat)
## [1] 0.04064889

What does each observation in the sampling distribution represent? ### Answer Each observation represents the sample proportion that makes up the sampling distribution.

How does the mean, standard error, and shape of the sampling distribution change as the sample size increases? ### Answer The standard error decreases as the sample size increases. This means as the sample size increases there is less variability. The shape of the sampling distribution appears normal as the sample increases.

How (if at all) do these values change if you increase the number of simulations? (You do not need to include plots in your answer.)

Answer

The values will change if the number of simulations is increased.

Exercise 7

Take a sample of size 15 from the population and calculate the proportion of people in this sample who think the work scientists do enhances their lives.

Using this sample, what is your best point estimate of the population proportion of people who think the work scientists do enchances their lives?

Answer

The best point estimate of the population proportion of people who think the work scientists do enchances their lives is 80%

sampl_p15 <- global_monitor %>%
                    rep_sample_n(size = 15, reps = 25, replace = TRUE) %>%
                    count(scientist_work) %>%
                    mutate(p_hat = n /sum(n)) %>%
                    filter(scientist_work == "Benefits")

sampl_p15
## # A tibble: 25 × 4
## # Groups:   replicate [25]
##    replicate scientist_work     n p_hat
##        <int> <chr>          <int> <dbl>
##  1         1 Benefits          13 0.867
##  2         2 Benefits          12 0.8  
##  3         3 Benefits          13 0.867
##  4         4 Benefits          13 0.867
##  5         5 Benefits          12 0.8  
##  6         6 Benefits          10 0.667
##  7         7 Benefits          12 0.8  
##  8         8 Benefits          12 0.8  
##  9         9 Benefits          12 0.8  
## 10        10 Benefits          14 0.933
## # … with 15 more rows
mean(sampl_p15$p_hat)
## [1] 0.816
sd(sampl_p15$p_hat)
## [1] 0.08877771

Exercise 8

Since you have access to the population, simulate the sampling distribution of proportion of those who think the work scientists do enchances their lives for samples of size 15 by taking 2000 samples from the population of size 15 and computing 2000 sample proportions. Store these proportions in as sample_props15.

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")

sample_props15
## # A tibble: 2,000 × 4
## # Groups:   replicate [2,000]
##    replicate scientist_work     n p_hat
##        <int> <chr>          <int> <dbl>
##  1         1 Benefits          10 0.667
##  2         2 Benefits           9 0.6  
##  3         3 Benefits          14 0.933
##  4         4 Benefits          13 0.867
##  5         5 Benefits          13 0.867
##  6         6 Benefits          14 0.933
##  7         7 Benefits          10 0.667
##  8         8 Benefits          11 0.733
##  9         9 Benefits          12 0.8  
## 10        10 Benefits          13 0.867
## # … with 1,990 more rows

Answer

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?

ggplot(data = sample_props15, aes(x = p_hat)) +
  geom_histogram(binwidth = 0.03) +
  labs(
    x = "p_hat (Benefit)",
    title = "Sampling distribution",
    subtitle = "Sample size = 15, Number of samples = 2000"
  )

Finally, calculate and report the population proportion.

 mean(sample_props15$p_hat)
## [1] 0.8022

Exercise 9

Change your sample size from 15 to 150, then compute the sampling distribution using the same method as above, and store these proportions in a new object called sample_props150.

Answer

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")

sample_props150
## # A tibble: 2,000 × 4
## # Groups:   replicate [2,000]
##    replicate scientist_work     n p_hat
##        <int> <chr>          <int> <dbl>
##  1         1 Benefits         120 0.8  
##  2         2 Benefits         119 0.793
##  3         3 Benefits         128 0.853
##  4         4 Benefits         115 0.767
##  5         5 Benefits         124 0.827
##  6         6 Benefits         128 0.853
##  7         7 Benefits         125 0.833
##  8         8 Benefits         120 0.8  
##  9         9 Benefits         116 0.773
## 10        10 Benefits         117 0.78 
## # … with 1,990 more rows
mean(sample_props150$p_hat)
## [1] 0.8016933

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

Based on this proportion, the mean is 80%.

ggplot(data = sample_props150, aes(x = p_hat)) +
  geom_histogram(binwidth = 0.02) +
  labs(
    x = "p_hat (Benefit)",
    title = "Sampling distribution",
    subtitle = "Sample size = 150, Number of samples = 2000"
  )

Exercise 10

Of the sampling distributions from 2 and 3, which has a smaller spread? If you’re concerned with making estimates that are more often close to the true value, would you prefer a sampling distribution with a large or small spread?

Answer

Comparing sampling distributions from 2 and 3 the spread with sample size 150 is smaller than the sampling distribution with sample size 15. I would prefer sampling a distribution with a smaller spread because there is less variability in the data. The smaller spread will give a closer approximation because the standard error will be smaller.