Here are my answers to the Week 4 Lab Activity of the couse Inferential Statistics with R presented by Coursera and conducted by Mine Çetinkaya-Rundel.

R Markdown

This is a R Markdown file. To a better viewing, it could be forked and knitted on RStudio to a html file or it could be viewed directly as a RPubs publication.

Packages

We will use the devtools package to install the statsr package associated with this course. We need to install and load this package.

install.packages("devtools")
library(devtools)

Now we can install the rest of the packages we will use during the course. Type the following commands in the Console as well:

install.packages("dplyr")
install.packages("ggplot2")
install.packages("shiny")
install_github("StatsWithR/statsr")

Getting Started

In August of 2012, news outlets ranging from the Washington Post to the Huffington Post ran a story about the rise of atheism in America. The source for the story was a poll that asked people, “Irrespective of whether you attend a place of worship or not, would you say you are a religious person, not a religious person or a convinced atheist?” This type of question, which asks people to classify themselves in one way or another, is common in polling and generates categorical data. In this lab we take a look at the atheism survey and explore what’s at play when making inference about population proportions using categorical data.

Load packages

In this lab we will explore the data using the dplyr package and visualize it using the ggplot2 package for data visualization. The data can be found in the companion package for this course, statsr.

Let’s load the packages.

library(statsr)
library(dplyr)
library(ggplot2)

The survey

The press release for the poll, conducted by WIN-Gallup International, can be accessed here.

Take a moment to review the report then address the following questions.

  1. How many people were interviewed for this survey?
    1. A poll conducted by WIN-Gallup International surveyed 51,000 people from 57 countries.
    2. A poll conducted by WIN-Gallup International surveyed 52,000 people from 57 countries.
    3. A poll conducted by WIN-Gallup International surveyed 51,917 people from 57 countries.
    4. A poll conducted by WIN-Gallup International surveyed 51,927 people from 57 countries.
      5. This poll surveyed 51,927 people from 57 countries.
  2. Which of the following methods were used to gather information?
    1. Face to face
    2. Telephone
    3. Internet
    4. All of the above
      5. A total of 51,927 persons were interviewed globally. In each country a national probability sample of around 1000 men and women was interviewed either face to face (35 countries; n=33,890), via telephone (11 countries; n=7,661) or online (11 countries; n=10,376).
  3. True / False: In the first paragraph, several key findings are reported. These percentages appear to be sample statistics.
    1. True
    2. False
      3. True. These percentages appear to be sample statistics derived from all the interviews conducted accross the 57 countries.
  4. True / False:The title of the report is “Global Index of Religiosity and Atheism”. To generalize the report’s findings to the global human population, We must assume that the sample was a random sample from the entire population in order to be able to generalize the results to the global human population. This seems to be a reasonable assumption.
    1. True
    2. False
      3. False. People were interviewed face to face, by fone and online. Strictly speaking, we could not assume a random sample since fone and online interviews are not necessarily independent.

The data

Turn your attention to Table 6 (pages 15 and 16), which reports the sample size and response percentages for all 57 countries. While this is a useful format to summarize the data, we will base our analysis on the original data set of individual responses to the survey. Load this data set into R with the following command.

data(atheism)
  1. What does each row of Table 6 correspond to?
    1. Countries
    2. Individual Persons
    3. Religions
      4. Each row of Table 6 corresponds to countries.
  2. What does each row of atheism correspond to?
    1. Countries
    2. Individual Persons
    3. Religions
      4. Each row of atheism correspond to individual persons. 
summary(atheism)
##              nationality           response          year     
##  Pakistan          : 5409   atheist    : 5498   Min.   :2005  
##  France            : 3359   non-atheist:82534   1st Qu.:2005  
##  Korea, Rep (South): 3047                       Median :2012  
##  Ghana             : 2995                       Mean   :2009  
##  Macedonia         : 2418                       3rd Qu.:2012  
##  Peru              : 2414                       Max.   :2012  
##  (Other)           :68390
head(atheism)
## # A tibble: 6 × 3
##   nationality    response  year
##        <fctr>      <fctr> <int>
## 1 Afghanistan non-atheist  2012
## 2 Afghanistan non-atheist  2012
## 3 Afghanistan non-atheist  2012
## 4 Afghanistan non-atheist  2012
## 5 Afghanistan non-atheist  2012
## 6 Afghanistan non-atheist  2012
tail(atheism)
## # A tibble: 6 × 3
##   nationality    response  year
##        <fctr>      <fctr> <int>
## 1     Vietnam non-atheist  2005
## 2     Vietnam non-atheist  2005
## 3     Vietnam non-atheist  2005
## 4     Vietnam non-atheist  2005
## 5     Vietnam non-atheist  2005
## 6     Vietnam non-atheist  2005

To investigate the link between these two ways of organizing this data, take a look at the estimated proportion of atheists in the United States. Towards the bottom of Table 6, we see that this is 5%. We should be able to come to the same number using the atheism data.

Create a new dataframe called us12 that contains only the rows in atheism associated with respondents to the 2012 survey from the United States:

us12 <- atheism %>%
  filter(nationality == "United States" , atheism$year == "2012")
  1. Next, calculate the proportion of atheist responses in the United States in 2012, i.e. in us12. True / False: This percentage agrees with the percentage in Table 6.
    1. True
    2. False
      3. True. Yes, it agrees with the Table 6. It is a sample statistic. 
us12_tab <- table(us12$response)
us12_tab
## 
##     atheist non-atheist 
##          50         952
p_hat <- us12_tab[[1]]/us12_tab[[2]]
p_hat
## [1] 0.05252101

Inference on proportions

As was hinted earlier, Table 6 provides sample statistics, that is, calculations made from the sample of 51,927 people. What we’d like, though, is insight into the population parameters. You answer the question, “What proportion of people in your sample reported being atheists?” with a statistic; while the question “What proportion of people on earth would report being atheists” is answered with an estimate of the parameter.

The inferential tools for estimating population proportion are analogous to those used for means in the last lab: the confidence interval and the hypothesis test.

Exercise: Write out the conditions for inference to construct a 95% confidence interval for the proportion of atheists in the United States in 2012. Are you confident all conditions are met?

The sampling distribution for \(\hat{p}\), taken from a sample of size \(n\) from a population with a true proportion \(p\), is nearly normal, i.e., allows inference when

  • the sample observations are independent and
  • we expected to see at least 10 successes and 10 failures in our sample, i.e. \(np \ge 10\) and \(n(1 − p)\ge 10\). This is the success-failure condition.

For confidence intervals, when we do no have information about the population parameter \(p\), we use the sample proportion \(\hat{p}\) to check the success-failure condition.

We have discussed the independence condition on question 4. We would need more information to state the independence condition for the United States survey, because we do not know if all the interviews in US were made face by face or if there were some made by telephone or internet.

On the other hand, the success-failure condition is clearly satisfied.

n_us = nrow(us12)
p_hat*n_us
## [1] 52.62605
(1-p_hat)*n_us
## [1] 949.3739
inference(y = response, data = us12, statistic = "proportion", type = "ci", method = "theoretical", success = "atheist")
## Single categorical variable, success: atheist
## n = 1002, p-hat = 0.0499
## 95% CI: (0.0364 , 0.0634)

Note that since the goal is to construct an interval estimate for a proportion, it’s necessary to specify what constitutes a `success'', which here is a response ofatheist`.

Although formal confidence intervals and hypothesis tests don’t show up in the report, suggestions of inference appear at the bottom of page 7: “In general, the error margin for surveys of this kind is \(\pm\) 3-5% at 95% confidence.”

Exercise: Imagine that, after reading a front page story about the latest public opinion poll, a family member asks you, “What is a margin of error?” In one sentence, and ignoring the mechanics behind the calculation, how would you respond in a way that conveys the general concept?

The margin of error quantifies the uncertainty about the result. The actual value is between the estimate minus the margin of error and the estimate plus the margin of error.

  1. Based on the R output, what is the margin of error for the estimate of the proportion of the proportion of atheists in US in 2012?
    1. The margin of error for the estimate of the proportion of atheists in the US in 2012 is 0.05.
    2. The margin of error for the estimate of the proportion of atheists in the US in 2012 is 0.025.
    3. The margin of error for the estimate of the proportion of atheists in the US in 2012 is 0.0135. 
n <- nrow(us12)
SE <- sqrt( p_hat *(1 - p_hat)/n  )
z_star <- qnorm(0.025, lower.tail = FALSE)
z_star*SE
## [1] 0.01381229

Note there is a slightly diffrence between the p_hat we have calculated and the p-hat given by the inference function.

Exercise: Using the inference function, calculate confidence intervals for the proportion of atheists in 2012 in two other countries of your choice, and report the associated margins of error. Be sure to note whether the conditions for inference are met. It may be helpful to create new data sets for each of the two countries first, and then use these data sets in the inference function to construct the confidence intervals.

br12 <- atheism %>%
  filter(nationality == "Brazil" , atheism$year == "2012")

inference(y = response, data = br12, statistic = "proportion", type = "ci", method = "theoretical", success = "atheist")
## Single categorical variable, success: atheist
## n = 2002, p-hat = 0.01
## 95% CI: (0.0056 , 0.0143)

ar12 <- atheism %>%
  filter(nationality == "Argentina" , atheism$year == "2012")

inference(y = response, data = ar12, statistic = "proportion", type = "ci", method = "theoretical", success = "atheist")
## Single categorical variable, success: atheist
## n = 991, p-hat = 0.0706
## 95% CI: (0.0547 , 0.0866)

How does the proportion affect the margin of error?

Imagine you’ve set out to survey 1000 people on two questions: are you female? and are you left-handed? Since both of these sample proportions were calculated from the same sample size, they should have the same margin of error, right? Wrong! While the margin of error does change with sample size, it is also affected by the proportion.

Think back to the formula for the standard error: \(SE = \sqrt{p(1-p)/n}\). This is then used in the formula for the margin of error for a 95% confidence interval: \(ME = 1.96\times SE = 1.96\times\sqrt{p(1-p)/n}\). Since the population proportion \(p\) is in this \(ME\) formula, it should make sense that the margin of error is in some way dependent on the population proportion. We can visualize this relationship by creating a plot of \(ME\) vs. \(p\).

The first step is to make a vector p that is a sequence from \(0\) to \(1\) with each number separated by \(0.01\). We can then create a vector of the margin of error (me) associated with each of these values of p using the familiar approximate formula (\(ME = 1.96 \times SE\)). Lastly, we plot the two vectors against each other to reveal their relationship.

d <- data.frame(p <- seq(0, 1, 0.01))
n <- 1000
d <- d %>%
  mutate(me = 1.96*sqrt(p*(1 - p)/n))
ggplot(d, aes(x = p, y = me)) +
  geom_line()

  1. Which of the following is false about the relationship between \(p\) and \(ME\).
    1. The \(ME\) reaches a minimum at \(p = 0\).
    2. The \(ME\) reaches a minimum at \(p = 1\).
    3. The \(ME\) is maximized when \(p = 0.5\).
    4. The most conservative estimate when calculating a confidence interval occurs when \(p\) is set to 1.
      5. The most conservative estimate when calculating a confidence interval occurs when \(p=0.5\).

The question of atheism was asked by WIN-Gallup International in a similar survey that was conducted in 2005. We assume here that sample sizes have remained the same. Table 4 on page 13 of the report summarizes survey results from 2005 and 2012 for 39 countries.

Answer the following two questions using the inference function. As always, write out the hypotheses for any tests you conduct and outline the status of the conditions for inference.

  1. True / False: There is convincing evidence that Spain has seen a change in its atheism index between 2005 and 2012.

    Hint: Create a new data set for respondents from Spain. Then use their responses as the first input on the inference, and use year as the grouping variable.
    1. True
    2. False
      3. To answer the question we set the following hypothesis test: \(H_0: p_{2005} = p_{2012}\); \(H_A : p_{2005} \ne p_{2012}\). 
spain <- atheism %>%
  filter(nationality == "Spain")

inference(y = response, x = as.factor(year), data = spain, statistic = "proportion", type = "ht", null = 0, success = "atheist", alternative = "twosided", method = "theoretical")
## Response variable: categorical (2 levels, success: atheist)
## Explanatory variable: categorical (2 levels) 
## n_2005 = 1146, p_hat_2005 = 0.1003
## n_2012 = 1145, p_hat_2012 = 0.09
## H0: p_2005 =  p_2012
## HA: p_2005 != p_2012
## z = 0.8476
## p_value = 0.3966

As the p-value is greater than the survey confidence level (\(\alpha > 0.05\)), we fail to reject the null hypothesis. In other words, the data does not provide enought evidence that Spain had seen a change in its atheism index between 2005 and 2012.

  1. True / False: There is convincing evidence that the United States has seen a change in its atheism index between 2005 and 2012.
    1. True
    2. False
      3. We set the same hypothesis test we used in the previous question, but now we filter the United States’ data. 
us <- atheism %>%
  filter(nationality == "United States")

inference(y = response, x = as.factor(year), data = us, statistic = "proportion", type = "ht", null = 0, success = "atheist", alternative = "twosided", method = "theoretical")
## Response variable: categorical (2 levels, success: atheist)
## Explanatory variable: categorical (2 levels) 
## n_2005 = 1002, p_hat_2005 = 0.01
## n_2012 = 1002, p_hat_2012 = 0.0499
## H0: p_2005 =  p_2012
## HA: p_2005 != p_2012
## z = -5.2431
## p_value = < 0.0001

As the p-value is lower than the survey confidence level (\(\alpha > 0.05\)), we reject the null hypothesis. In other words, the data provides enought evidence that the United States had seen a change in its atheism index between 2005 and 2012.

  1. If in fact there has been no change in the atheism index in the countries listed in Table 4, in how many of those countries would you expect to detect a change (at a significance level of 0.05) simply by chance?

    Hint: Type 1 error.
    1. 0
    2. 1
    3. 1.95
    4. 5
      5. A Type 1 error is rejecting the null hypothesis when \(H_0\) is actually true. As a general rule of thumb, for those cases where the null hypothesis is actually true, we do not want to incorrectly reject \(H_0\) more than 5% of the time, which is in fact the case in this context, because the significance level is 0.05. Since there are 39 countries, we would expect to detect a change simply by chance in approximately 39 \(\times\) 0.05 countries. 
39*0.05
## [1] 1.95
  1. 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. Do not use the data set to answer this question.
    1. 2401 people
    2. At least 2401 people
    3. 9604 people
    4. At least 9604 people
      5. We use the expression: \[ ME = z^{*}\sqrt{ \frac{p(1-p)}{n}}\] 
p_maxME = 0.5
z_star =  qnorm(0.025, lower.tail = FALSE)
ME = 0.01
n_min = z_star*z_star*p_maxME*(1-p_maxME)/(ME*ME)
n_min
## [1] 9603.647

At least 9604 people.

This is a product of OpenIntro that is released under a Creative Commons Attribution-ShareAlike 3.0 Unported. This lab was written for OpenIntro by Andrew Bray and Mine Çetinkaya-Rundel.