Chapter 6 Inference for Categorical Data

load the package

library(DATA606) 
## 
## Welcome to CUNY DATA606 Statistics and Probability for Data Analytics 
## This package is designed to support this course. The text book used 
## is OpenIntro Statistics, 3rd Edition. You can read this by typing 
## vignette('os3') or visit www.OpenIntro.org. 
##  
## The getLabs() function will return a list of the labs available. 
##  
## The demo(package='DATA606') will list the demos that are available.
## 
## Attaching package: 'DATA606'
## The following object is masked from 'package:utils':
## 
##     demo
library(ggplot2)
## Warning: package 'ggplot2' was built under R version 3.3.3
library(knitr)

6.7.1 Inference for a single proportion

6.6 2010 Healthcare Law. On June 28, 2012 the U.S. Supreme Court upheld the much debated 2010 healthcare law, declaring it constitutional. A Gallup poll released the day after this decision indicates that 46% of 1,012 Americans agree with this decision. At a 95% confidence level, this sample has a 3% margin of error. Based on this information, determine if the following statements are true or false, and explain your reasoning.

  1. We are 95% confident that between 43% and 49% of Americans in this sample support the decision of the U.S. Supreme Court on the 2010 healthcare law.
# False. confidence interval only attempt to capture population parameters. We are 95% confident that the fraction of the population (but not the sample) supporting the decision is between 43% and 49%. 
  1. We are 95% confident that between 43% and 49% of Americans support the decision of the U.S. Supreme Court on the 2010 healthcare law.
# True. This is the definition of confidence interval.
  1. If we considered many random samples of 1,012 Americans, and we calculated the sample proportions of those who support the decision of the U.S. Supreme Court, 95% of those sample proportions will be between 43% and 49%.
# True. This is the definition of confidence interval.
  1. The margin of error at a 90% confidence level would be higher than 3%.
# False. Margion of error = Z_star * Standard Error. Z score for 90% power is 1.65 which is less than Z score for 95% power, 1.96. So the margion of error will become samller than 3%.

6.12 Legalization of marijuana, Part I. The 2010 General Social Survey asked 1,259 US residents: “Do you think the use of marijuana should be made legal, or not?” 48% of the respondents said it should be made legal.

  1. Is 48% a sample statistic or a population parameter? Explain.
# 48% is a sample statistic becase this number is derived from a sample of 1259 US redidents but not from the US population.
  1. Construct a 95% confidence interval for the proportion of US residents who think marijuana should be made legal, and interpret it in the context of the data.
# Is distribution of p being nearly normal?
# 1. the sample observations are independent.
# 2. 48% of the respondents said it should be made legal. So at least 10 successes and 10 failures in our sample.
# So The tribution of sample proportion p being nearly normal
# For 95% confidence interval, alpha = 0.05, Z-score = 1.96
Z1 <- 1.96
p <- 0.48
n1 <- 1259
margion_error1 <- Z1 * sqrt(p*(1-p)/n1)
ci_min <- (p - margion_error1)
ci_max <- (p + margion_error1)
paste("The 95% confidence interval for the proportion of US residents who think marijuana should be made legal range from", round(ci_min*100,2), "%", "to", round(ci_max*100,2), "%.", "We are 95% confident that propotion of US population that support the legal use of marijuana range from", round(ci_min*100,2), "%", "to", round(ci_max*100,2), "%.")
## [1] "The 95% confidence interval for the proportion of US residents who think marijuana should be made legal range from 45.24 % to 50.76 %. We are 95% confident that propotion of US population that support the legal use of marijuana range from 45.24 % to 50.76 %."
  1. A critic points out that this 95% confidence interval is only accurate if the statistic follows a normal distribution, or if the normal model is a good approximation. Is this true for these data? Explain.
suc <- 1259 * p
fail <- 1259 * (1-p)
print(suc) 
## [1] 604.32
print(fail) 
## [1] 654.68
# 1. Observations are independent. The poll is based on a simple random sample and consists of fewer than 10% of the US residents, which verifies independence.
# 2. Success-failure condition. The sample size must is fairly large, which is checked using the success-failure condition. There were 1259 * p = 604 "successes" and 1259 * (1-p) = 654 "failures" in the sample, both easily greater than 10.

# So the sampling distribution of p is nearly normal.
  1. A news piece on this survey’s findings states, “Majority of Americans think marijuana should be legalized.” Based on your confidence interval, is this news piece’s statement justified?
# False. The confidence interval is 45.24 % to 50.76 %, which barely reach 50%. So it is not correct to say "Majority of Americans think marijuana should be legalized."

6.20 Legalize Marijuana, Part II. As discussed in Exercise 6.12, the 2010 General Social Survey reported a sample where about 48% of US residents thought marijuana should be made legal. If we wanted to limit the margin of error of a 95% confidence interval to 2%, about how many Americans would we need to survey?

p3 <- 0.48
margion_error3 <- 0.02
# For 95% power, Z score = 1.96
Z3 <- 1.96
# Because Margion Error = Z_star * SE
SE3 <- margion_error3/Z3
# Because SE = sqrt(p*(1-p)/n)
n3 <- p3*(1-p3)/SE3^2

paste("If we wanted to limit the margin of error of a 95% confidence interval to 2%, we will need to survay", round(n3+0.5),"Americans.")
## [1] "If we wanted to limit the margin of error of a 95% confidence interval to 2%, we will need to survay 2398 Americans."

6.7.2 Difference of two proportions

6.28 Sleep deprivation, CA vs. OR, Part I. According to a report on sleep deprivation by the Centers for Disease Control and Prevention, the proportion of California residents who reported insufficient rest or sleep during each of the preceding 30 days is 8.0%, while this proportion is 8.8% for Oregon residents. These data are based on simple random samples of 11,545 California and 4,691 Oregon residents. Calculate a 95% confidence interval for the difference between the proportions of Californians and Oregonians who are sleep deprived and interpret it in context of the data.

# NULL hypothesis: There is no difference in proportions of Californians and Oregonians who are sleep deprived. mu_CA = mu_OR 
# Alternative hypothesis: There is difference in proportions of Californians and Oregonians who are sleep deprived. mu_CA != mu_OR 
n_CA <- 11545
p_CA <- .08
n_OR <- 4691
p_OR <- .088
# For 95% power, Z score = 1.96
Z4 <- 1.96

SE_CA <- sqrt(p_CA *(1-p_CA)/n_CA)
ME_CA  <- Z4*SE_CA
ci_min_CA <- p_CA - ME_CA
ci_max_CA <- p_CA + ME_CA

SE_OR <- sqrt(p_OR *(1-p_OR)/n_OR)
ME_OR <- Z4*SE_OR
ci_min_OR <- p_OR - ME_OR
ci_max_OR <- p_OR + ME_OR
paste("95% confidence interval for proportion of Sleep deprivation in California is", round(ci_min_CA*100,2), "%", "to", round(ci_max_CA*100,2), "%.")
## [1] "95% confidence interval for proportion of Sleep deprivation in California is 7.51 % to 8.49 %."
paste("95% confidence interval for proportion of Sleep deprivation in Oregon is", round(ci_min_OR*100,2), "%", "to", round(ci_max_OR*100,2), "%.")
## [1] "95% confidence interval for proportion of Sleep deprivation in Oregon is 7.99 % to 9.61 %."
SE_dif <- sqrt(p_CA *(1-p_CA)/n_CA + p_OR *(1-p_OR)/n_OR)
mu_dif <- p_CA - p_OR
# calculate Z-score
Z5 <- (mu_dif-0)/SE_dif
pVal5 <- 2*(1-pnorm(abs(Z5)))

if(pVal5 > 0.05){
  paste("p-Value =", round(pVal5,4), "> 0.05",  ", there is no strong evidence to reject H0.")
}else{
  paste("p-Value =", round(pVal5,4), "< 0.05", ", so we reject H0. There is difference in proportions of Californians and Oregonians who are sleep deprived.")
}
## [1] "p-Value = 0.0988 > 0.05 , there is no strong evidence to reject H0."

** 6.7.3 Testing for goodness of fit using chi-square**

6.44 Barking deer. Microhabitat factors associated with forage and bed sites of barking deer in Hainan Island, China were examined from 2001 to 2002. In this region woods make up 4.8% of the land, cultivated grass plot makes up 14.7% and deciduous forests makes up 39.6%. Of the 426 sites where the deer forage, 4 were categorized as woods, 16 as cultivated grassplot, and 61 as deciduous forests. The table below summarizes these data.

df_deer <- data.frame(Woods = 4, Cultivated_grassplot = 16, Deciduous_forests = 67, Other = 345, Total = 426)
df_deer
##   Woods Cultivated_grassplot Deciduous_forests Other Total
## 1     4                   16                67   345   426
  1. Write the hypotheses for testing if barking deer prefer to forage in certain habitats over others.
# NULL hypothesis: Barking deer do not prefer to forage in certain habitats over others. 
# Alternative hypothesis: Barking deer do prefer to forage in certain habitats over others.  
  1. What type of test can we use to answer this research question?
# We can use a chi square test.
  1. Check if the assumptions and conditions required for this test are satisfied.
# There are 4 cases inthe woods sectionn, but the expected cases of Woods is as following:
exp_W <- 426 * .048
exp_W
## [1] 20.448
# Assumming that the deer are not influenced by each other and they are independent of all the other cases in the table. 
# Each particular scenario have at least 5 expected cases.
# So the two conditions is satisfied for performing a chi-square test.
  1. Do these data provide convincing evidence that barking deer prefer to forage in certain habitats over others? Conduct an appropriate hypothesis test to answer this research question.
exp_C <- 426 * .147
exp_D <- 426 * .396
exp_O <- 426 * (1-.048-.147-.396)

chi_deer <- (4-exp_W)^2/exp_W + (16-exp_C)^2/exp_C + (67-exp_D)^2/exp_D + (345-exp_O)^2/exp_O 

paste("The Chi-square for deer preference for certain habitats is", round(chi_deer,2))
## [1] "The Chi-square for deer preference for certain habitats is 276.61"
pVal7 <- dchisq(chi_deer, df = 4-1)
if(pVal7 > 0.05){
  paste("p-Value =", pVal7, "> 0.05",  ", there is no strong evidence to reject H0.")
}else{
  paste("p-Value =", pVal7, "< 0.05", ", so we reject H0. Barking deer do prefer to forage in certain habitats over others.")
}
## [1] "p-Value = 5.70144101684255e-60 < 0.05 , so we reject H0. Barking deer do prefer to forage in certain habitats over others."

6.7.4 Testing for independence in two-way tables

6.48 Coffee and Depression. Researchers conducted a study investigating the relationship between caffeinated coffee consumption and risk of depression in women. They collected data on 50,739 women free of depression symptoms at the start of the study in the year 1996, and these women were followed through 2006. The researchers used questionnaires to collect data on ca???einated coffee consumption, asked each individual about physician-diagnosed depression, and also asked about the use of antidepressants. The table below shows the distribution of incidences of depression by amount of caffeinated coffee consumption.

a8 <- c(670,11545,12215)
b8 <- c(373,6244,6617) 
c8 <- c(905,16329,17234)
d8 <- c(564,11726,12290)
e8 <- c(95,2288,2383)
f8 <- c(2607,48132,50739)
df_coffee <- data.frame(cup_1_week =a8 , cup_2_6_week = b8, cup_1_day = c8, cup_2_3_day = d8, cup_4_day = e8,   Total = f8)
row.names(df_coffee) <- c("Clinical_depresssion_Yes", "Clinical_depresssion_No", "Clinical_depresssion_Total")
kable(df_coffee)
cup_1_week cup_2_6_week cup_1_day cup_2_3_day cup_4_day Total
Clinical_depresssion_Yes 670 373 905 564 95 2607
Clinical_depresssion_No 11545 6244 16329 11726 2288 48132
Clinical_depresssion_Total 12215 6617 17234 12290 2383 50739
  1. What type of test is appropriate for evaluating if there is an association between coffee intake and depression?
# Two table chi-square test is the appropriate test.
  1. Write the hypotheses for the test you identified in part (a).
# NULL hypothesis: the proportion of women experienced clinical depression are the same no matter how much coffee they drink.
# Alternative hypothesis: the proportion of women experienced clinical depression are different among women who drink diffrent amount of coffee.
  1. Calculate the overall proportion of women who do and do not suffer from depression.
p_dep <- df_coffee[1,6]/df_coffee[3,6]
paste("The proportion of women who suffer from depression is", round(100*p_dep,2),"%.")
## [1] "The proportion of women who suffer from depression is 5.14 %."
p_notdep <- df_coffee[2,6]/df_coffee[3,6]
paste("The proportion of women who do not suffer from depression is", round(100*p_notdep,2),"%.")
## [1] "The proportion of women who do not suffer from depression is 94.86 %."
  1. Identify the expected count for the highlighted cell, and calculate the contribution of this cell to the test statistic, i.e. (Observed ??? Expected)^2/Expected.
# 373 is highlighted
exp_2_6_pw <-  df_coffee[3,2] * p_dep
paste("Expected value of 2 to 6 cups per week and depressed women is:", exp_2_6_pw)
## [1] "Expected value of 2 to 6 cups per week and depressed women is: 339.985395849347"
test_stat_26 <- (373 -exp_2_6_pw)^2/exp_2_6_pw
paste("the test statistic of 2 to 6 cups per week and depressed women is:", test_stat_26)
## [1] "the test statistic of 2 to 6 cups per week and depressed women is: 3.20591443200495"
  1. The test statistic is #2 = 20.93. What is the p-value?
# degree of freedom is 
df9 <- (5-1) * (2-1)
pVal9 <- dchisq(20.93, df = df9)
paste("p-Value:", pVal9)
## [1] "p-Value: 0.000149216718130047"
  1. What is the conclusion of the hypothesis test?
if(pVal9 > 0.05){
  paste("p-Value =", pVal9, "> 0.05",  ", there is no strong evidence to reject H0.")
}else{
  paste("p-Value =", pVal9, "< 0.05", ", so we reject H0. The proportion of women experienced clinical depression are different among women who drink diffrent amount of coffee.")
}
## [1] "p-Value = 0.000149216718130047 < 0.05 , so we reject H0. The proportion of women experienced clinical depression are different among women who drink diffrent amount of coffee."
  1. One of the authors of this study was quoted on the NYTimes as saying it was “too early to recommend that women load up on extra coffee” based on just this study. Do you agree with this statement? Explain your reasoning.
# I agree that it is "too early to recommend that women load up on extra coffee" based on just this study. Because this is a observation study, the results do not necessitate a causation. Besides, there are many other confounding variables could influence the results.