BIOS:4120 - Lab 11: Fisher Test, CI for OR, and Quiz 3 Review
April 2-3, 2024
Objectives
- Revisit Fisher’s Exact Test
- Discuss the confidence intervals for the odds ratio
- Quiz 3 Review
Note: The material in the first part of this lab will not be on quiz 3.
Fisher Test and Confidence Interval
Contingency tables with Lister data
Let’s start by reviewing the lister dataset.
lister <- read.delim('https://github.com/IowaBiostat/data-sets/raw/main/lister/lister.txt')
table(lister) # make table of the data## Outcome
## Group Died Survived
## Control 16 19
## Sterile 6 34Do not fret that this table seems like it’s in reverse. It turns out that the table being in this “reverse” order doesn’t actually matter for when we do our analysis, but to look at the data as it was displayed during lecture, we can switch the row and column order with a bit of manipulation to the dataset.
lister$Group <- factor(lister$Group, levels = c("Sterile", "Control")) # specify order you want in the levels
lister$Outcome <- factor(lister$Outcome, levels = c("Survived", "Died"))
(lister_tab <- table(lister)) # we will use this table in later analysis## Outcome
## Group Survived Died
## Sterile 34 6
## Control 19 16Alternatively, if we knew the counts of the data themselves, we could make a table in the following way:
lister_manual <- rbind(c(34, 6),
c(19, 16))
# Add row and column names
rownames(lister_manual) <- c("Sterile", "Control")
colnames(lister_manual) <- c("Survived", "Died")
lister_manual## Survived Died
## Sterile 34 6
## Control 19 16# Alternative way to do this:
lister_manual2 <- data.frame("Group" = c(rep("Control", 35),
rep("Sterile", 40)),
"Outcome" = c(rep("Survived", 19),
rep("Died", 16),
rep("Survived", 34),
rep("Died", 6)))
# specify order of rows and columns in the levels:
lister_manual2$Group <- factor(lister_manual2$Group, levels = c("Sterile", "Control"))
lister_manual2$Outcome <- factor(lister_manual2$Outcome, levels = c("Survived", "Died"))
table(lister_manual2)## Outcome
## Group Survived Died
## Sterile 34 6
## Control 19 16Note: You may need to do something like this for your next homework assignment
Fisher’s exact test and odds ratio
Remember that Fisher’s Exact Test can be performed using the fisher.test() function. In lab last week using the lister dataset, we saw that this function not only gives you the p-value but also the odds ratio (OR) and a 95% CI for the odds ratio.
fisher.test(lister_tab)##
## Fisher's Exact Test for Count Data
##
## data: lister_tab
## p-value = 0.005018
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
## 1.437621 17.166416
## sample estimates:
## odds ratio
## 4.666849Note that this Fisher’s Exact Test odds ratio is calculated using the formula:
\(\frac{a/b}{c/d}=\frac{ad}{bc}\) , where a, b, c, and d are defined by their location according to the following table:
## Success Failure
## Option 1 a b
## Option 2 c dThe easiest way to think of this odds ratio is to say:
If OR > 1: The odds of success are (OR-1)*100 percent higher if we do Option 1 than if we do Option 2.
If OR < 1: The odds of success are (1-OR)*100 percent lower if we do Option 1 than if we do Option 2.
If OR = 1: The odds of success are equal between Option 1 and Option 2.
Be careful when calculating the OR that you interpret what you calculated. For instance, if I were to switch my Success and Failure column names, I’d be calculating the OR for Failure given Option 1 instead of Success. If you switch both the rows and column simultaneously (like we did with the very first table), the odds ratio will be unchanged.
Confidence interval for odds ratio
Calculating a confidence interval for an odds ratio is slightly more complicated than what we’ve done so far because it’s on a different scale. Here are the steps:
Find the odds ratio: \(\frac{ad}{bc}\)
Find the log odds ratio (natural log): \(\log{(OR)}\)
Find the error term: \(SE_{logOR}=\sqrt{\frac{1}{a} + \frac{1}{b} + \frac{1}{c} + \frac{1}{d}}\)
Calculate the CI on this scale: \(CI_{logOR}=\log{(OR)} \pm z_{\alpha/2} * SE_{logOR}\)
Calculate the CI on the original scale: \(CI = \exp{(CI_{logOR})}\)
Let’s find a CI for the OR for survival given sterile procedure by hand.
## Outcome
## Group Survived Died
## Sterile 34 6
## Control 19 16## Outcome
## Group Survived Died
## Sterile a b
## Control c d- Find the OR.
OR <- (34*16) / (19*6)
OR## [1] 4.77193- Find the log(OR)
logOR <- log(OR) # the log function in R is the natural log
logOR## [1] 1.562751- Find the error term.
SElogOR <- sqrt((1 / 19) + (1 / 6) + (1 / 16) + (1 / 34))
SElogOR## [1] 0.557862- Calculate the CI on this scale.
logCI <- logOR + qnorm(c(0.025, 0.975)) * SElogOR
logCI## [1] 0.4693614 2.6561402- Calculate the CI on the original scale.
CI <- exp(logCI)
CI## [1] 1.598973 14.241215Interpretation: We can say with 95% confidence that the true odds ratio for survival with sterile procedure relative to the control procedure is (1.599, 14.241), indicating significantly increased odds of survival on the sterile procedure compared to the non-sterile procedure. (Significant because the CI does not include the value 1.)
Note that this is the same as the relative odds of dying with the control surgery. This happens because the OR is symmetric.
Quiz 3 Review
A flowchart for when to use each distribution
General approach to stats problems
Step 1: Figure out and list the given information.
Write out the values of \(n, \bar{x}, s\), etc. for everything given in the prompt.
Step 2: Figure out what the problem is asking. Write the null hypothesis (if applicable).
What is the population we’re looking at? What kind of distributions can be used with the information given? Which is most appropriate in this situation? (see chart) Do we want a test or an interval? Most importantly, what is the question being asked?
If Defining Hypotheses: The null case is the one based on no change, difference, or improvement. The alternative is what we want to show and that there is a change or difference. Remember that the hypotheses describe and apply to a population parameter.
Step 3: Write out the equation for the test statistic or interval.
In this class, most test statistics will be similar to the form: \(z=\frac{(\hat{X}-X_0)}{SE_{\hat{X}}}\), and most confidence intervals will look similar to this: \(\hat{X} \pm z^**SE_{\hat{X}}\),
where:
- \(X\) is the parameter of
interest
- \(\hat{X}\) is the estimate of the
parameter
- \(X_0\) is the value of the
parameter under the null hypothesis
- \(SE_{\hat{X}}\) is the standard
error of the estimate
- \(z^*\) is the critical value for the interval (either from the z or t distribution)
Step 4: Plug in the values.
Calculate the statistic or interval using the given values. Find the p-value (if applicable). This means that if we’re conducting a test, compare the test statistic to the distribution we picked in step 2, and find the probability of being as or more extreme than the value we observed.
Step 5: Interpret in the context of the problem.
This might very well be the most important step. Our interpretation is dependent on our study, so it varies. Below are some examples.
Situation: Paired study looking at the effect of oat bran consumption compared to cornflakes levels where the calculated p-value was 0.005
- “There is strong evidence that eating oat bran as opposed to corn flakes leads to lower cholesterol levels.”
Situation: Weights of U.S. Males (mean believed to be 172.2) where observed mean was 180, and the calculated p-value was 0.09.
- “This study provides borderline evidence that the true mean weight of males in the United States is greater than 172.2 pounds.”
Table 1. Modifiers based on p-value:
| p | Evidence against null |
|---|---|
| > .1 | insufficient |
| .1 | borderline |
| .05 | moderate |
| .01 | strong |
| .001 | overwhelming |
Concepts to know
Effect size
Recall the Nexium/Prilosec example. If we have a large enough n, we can find statistical significance in the smallest true difference in means. Here, the difference, while real, was only 3%. That’s effectively REALLY close to being the same. For practical purposes, this difference doesn’t matter, but the p-value doesn’t tell us that the actual difference is little, just that it is significant.
Paired data
Paired data can be analyzed using the binomial distribution (proportion of patients that saw any improvement) or using continuous methods (using 1-sample methods on a set of the differences between the two).
Power
Given that the null hypothesis is false, the power of a test is the probability that we will correctly reject the null.
Properties of power:
- Power = (1 - \(\beta\)), where \(\beta\) is the type II error rate
- Power increases as sample size increases
- Power increases as type I error (\(\alpha\)) increases
- Power decreases as standard deviation increases
- Power increases as effect size increases
Practice Problems
1. We are interested in testing whether a certain at-risk population for diabetes has a daily sugar intake that is equal to the general population, which is equal to 77 grams/day. A sample of size 37 was taken from this at-risk population, and we obtained a sample mean of 80 and sample standard deviation of 11 grams. Perform a hypothesis test to test whether this population has a significantly different mean sugar intake from 77 grams.
Answer
We have continuous data and we have the sample standard deviation(s), so we should perform a t-test. The null hypothesis is that the true mean daily sugar intake is 77 grams/day.
\[ \bar{x} = 80\\ \mu_0 = 77\\ n = 37\\ s = 11\\ df = 37-1=36 \]
\[ SE = \frac{s}{\sqrt{n}} = \frac{11}{\sqrt{37}} = 1.808\\ t = \frac{\bar{x}}{SE} = \frac{80-77}{1.808} = 1.659 \]
The p-value is the probability of observing a value as extreme or more extreme than the one we observed under the null hypothesis. To compute it, find 2*Pr(T > 1.659). Using the t table, look at the row corresponding to our df. Since 36 is not on the table, you can just round to the nearest df you see on the table (or always round down if you like).
In that case, let’s look at df = 35, our test statistic t = 1.66 is between 1.47 and 1.69 which corresponds to the alpha values 0.15 and 0.1 so we can report the p-value in the following way:
0.1 < p-value < 0.15
Recall that the t table Patrick provides gives the area OUTSIDE the t-value so the p-value directly corresponds to the alpha values at the top of the table
We conclude that we do not have significant evidence to conclude that the at-risk mean sugar intake is different from the general population mean sugar intake of 77 grams/day (0.1 < p < 0.15).
mu <- 77
xbar <- 80
s <- 11
n <- 37
# Use a t-test for the sample
t <- (xbar-mu) / (s / sqrt(n))
2*pt(abs(t), n-1, lower.tail=FALSE)## [1] 0.10581772. The distribution of LDL cholesterol levels in a certain population is approximately normal with mean 90 mg/dl and standard deviation 8 mg/dl.
a) What is the probability an individual will have a LDL cholesterol level above 95 mg/dl?
Answer
The data are approximately normal and we are looking at INDIVIDUALS.
\[ x = 95 \\ \mu = 90 \\ \sigma = 8 \\ P(Z > \frac{x - \mu}{\sigma}) \\ = P(Z > \frac{95-90}{8}) = P(Z > 0.625) \\ = 1 - P(Z < 0.625) = P(Z <-0.625) \\ = 1 - 0.736 = 0.264 \]
mu <- 90
sigma <- 8
(p <- pnorm(95, mu, sigma, lower.tail=FALSE))## [1] 0.2659855b) Suppose we have a sample of 10 people from this population. What is the probability of exactly 3 of them being above 95 mg/dl?
Answer
We now can use the binomial distribution to find this probability where \(\pi\) is the probability of an INDIVIDUAL being above 95 mg/dl. We found exactly this probability in part a.
\[ X \sim Binom(n = 10, \pi = 0.264) \\ p(x) = \frac{n!}{x!(n-x)!}\pi^x(1-\pi)^{n-x} \\ p(3) = \frac{10!}{3!(7)!}0.264^3(1-0.264)^{7} \\ p(3) = 0.258 \]
dbinom(3, 10, p)## [1] 0.2592314c) Take the sample of size 10, as in part b. What is the probability that the sample mean will be above 95 mg/dl?
Answer
The data are still approximately normal and we are interested in the SAMPLE MEAN being above 95 mg/dl. This means we need to use the Z formula but use SE instead of the population sd (\(\sigma\))
\[ x = 95 \\ \mu = 90 \\ \sigma = 8 \\ n = 10 \\ SE = \frac{\sigma}{\sqrt{n}} = \frac{8}{\sqrt{10}} = 2.529822 \\ P(Z > \frac{x - \mu}{SE}) \\ = P(Z > \frac{95-90}{2.529822}) = P(Z > 1.976) \\ = 1 - P(Z < 1.976) = P(Z < -1.976) \\ = 1 - 0.976 = 0.024 \]
z <- (95-90)/(sigma/sqrt(10))
(p <- pnorm(z, lower.tail=FALSE))## [1] 0.02405341d) Suppose we take 5 samples of size 10 from the population. What is the probability that at least one of the sample means will be greater than 95 mg/dl?
Answer
Similarly to part b, we can use the Binomial distribution. Now the probability of succes (\(\pi\)) corresponds to the probability that a SAMPLE MEAN of sample size 10 will be greater than 95 mg/dl. Hopefully, you see that we calculated this exact probability in part c! Be careful here, the n for the Binomial distribution is how many samples we take, not the size of the samples.
\[ X \sim Binom(n = 5, \pi = 0.024) \\ p(x) = \frac{n!}{x!(n-x)!}\pi^x(1-\pi)^{n-x} \\ P(X \geq 1) = 1 - P(X \leq 0) = 1 - P(X = 0) = 1-p(0) \\ 1 - p(0) = 1- (\frac{5!}{0!(5)!}0.024^0(1-0.024)^{5}) \\ 1-p(0) = 1 - (1-0.024)^5 = 1 - 0.886 = 0.114 \]
Theh probability that AT LEAST ONE of the SAMPLE MEANS will be greater than 95 mg/dl is 0.114.
1-dbinom(0, 5, p)## [1] 0.11461893. In the city of Chicago, 235 baseball fans were sampled and asked if they cheer for the White Sox or the Cubs. 155 of the questioned people preferred the Cubs, while the remaining 80 fans preferred the White Sox. Use normal approximation to test the hypothesis that the Cubs and White Sox have an equal number of fans in the city. Construct a 95% Confidence Interval for the true proportion of Chicago fans who cheer for the Cubs, still using normal approximation.
Answer
We have binary/proportion data and we want to use a normal approximation. The null hypothesis is that each team has equal number of fans in the city, this corresponds to a proportion of 0.5.
\[ \pi_0 = 0.5 \\ n = 235 \\ \hat\pi = \frac{155}{235} = 0.6595745 \\ SE = \sqrt{\frac{\pi_0(1-\pi_0)}{n}} = \sqrt{\frac{0.5(1-0.5)}{235}} = 0.0326164 \]
\[ z = \frac{\hat\pi - \pi_0}{SE} = \frac{0.6595745 - 0.5} {0.0326164} = 4.89 \\ p = 2* P(Z > 4.89) = 2 * P(Z < -4.89) \approx 0.000 \]
# Writing down our known values
n <- 235
p <- 155 / n
p_0 <- 0.5
alpha <- 0.05
# Finding our z-statistic
z <- (p - p_0) / (sqrt(p_0*(1 - p_0) / n))
# Using our z-statistic to get our p-value
pvalue <- 2*(1 - pnorm(z, mean = 0, sd = 1))
pvalue## [1] 9.958308e-07There is overwhelming evidence to suggest that the proportion of Cubs fans is not equal to the proportion of White Sox fans in Chicago.
Now we use the corresponding confidence interval for proportions.
NOTE: SE for the confidence interval is different than the SE for the hypothesis test when working with the normal approximation for proportions!!
\[ n = 235 \\ \hat\pi = \frac{155}{235} = 0.6595745\\ SE = \sqrt{\frac{\hat\pi(1-\hat\pi)}{n}} = \sqrt{\frac{0.6596(1-0.6596)}{235}} = 0.03091015 \\ z_{\alpha/2} = z_{0.05/2} = z_{0.025} = -1.96 \]
95% CI
\[ \hat\pi \pm z_{\alpha/2} * SE \\ 0.6595745 \pm (1.96 * 0.03091015) \\ (0.599, 0.720) \]
SE <- sqrt(p*(1-p) / n)
z_alpha_2 <- qnorm(0.025)
p + c(1, -1)*z_alpha_2*SE## [1] 0.5989906 0.7201584Our null hypothesis value \(\pi_0=0.5\) is not contained in the interval. This is further confirmation that we can reject the null hypothesis, as it is not in the range of plausible values.
4. A pharmaceutical company is testing the effectiveness of a new drug for lowering blood pressure in people with prehypertension (systolic BP between 120 mmHg and 139 mmHg). They recruit 9 participants and measure their systolic blood pressure before and after taking the drug for one month. The mean systolic blood pressure for participants before taking the drug was 125 mmHg. The mean systolic blood pressure one month after taking the drug was 120.8 mmHg. The standard deviation for the difference (after - before) was 3.4 mmHg. This information is summarized in the table below.
| Before drug | After drug | Difference | |
|---|---|---|---|
| Mean | 125 | 120.8 | |
| Sd | 3.9 | 3.1 | 3.4 |
a) Test the hypothesis that the drug has no effect on lowering systolic blood pressure and state your conclusion.
Answer
The data are paired here because each participant’s blood pressure is measured both before and after taking the drug. Since we are given the means of the two groups and the sample standard deviation of the difference, we can use the t-test on the difference (use the formula just like we usually do but make \(\bar{x}\) the difference in means). The difference in means (after - before) is -4.2 mmHg. The null hypothesis is that the drug has no effect on lowering systolic blood pressure which is equivalent to \(\mu_0 = 0\)
\[ \bar{x} = -4.2 \\ \mu_0 = 0 \\ s = 3.4 \\ n = 9 \\ df = 8 \\ SE = \frac{s}{\sqrt{n}} = \frac{3.4}{\sqrt{9}} = 1.133\\ t = \frac{\bar{x} - \mu_0}{SE} = \frac{-4.2}{1.133} = -3.706 \]
To find the p-value, we can look for the +3.706 in the row with 8 degrees of freedom on the t table. This value falls between 3.36 and 3.83 which corresponds to a p-value between 0.01 and 0.005.
There is strong evidence to suggest that the difference in blood pressures before and after treatment is not equal to 0 (0.005 < p < 0.01). From the data we observed, we can conclude that the treatment is helpful in lowering blood pressure.
n <- 9
s <- 3.4
SE <- 3.4 / sqrt(n)
(t <- (-4.2 - 0) / SE)## [1] -3.705882(p_value <- 2*pt(t, 8))## [1] 0.005991706b) Calculate the 95% confidence interval for the mean difference in blood pressure before and after taking the drug.
Answer
We have all the parts we need except for \(t_{\alpha/2, n-1}\). To find this using the table, go to the df = 8 row and \(\alpha = 0.05\) column. We don’t look up \(\alpha = 0.025\) because the table Patrick provides accounts for both tails but the value we get still corresponds to the 0.025 quantile. When we do that, we get the corresponding t-value is 2.31.
\[ \bar{x} = -4.2 \\ s = 3.4 \\ n = 9 \\ df = 8 \\ SE = \frac{s}{\sqrt{n}} = \frac{3.4}{\sqrt{9}} = 1.133\\ t_{0.025, 8} = 2.31 \\ 95\% CI: \bar{x} \pm t_{0.025, 8}*SE \\ -4.2 \pm 2.31 * 1.133 \\ (-6.82, -1.58) \]
t_0.025 <- qt(0.975, 8)
-4.2 + c(-1,1)*t_0.025*SE## [1] -6.813471 -1.5865295. An article in the New England Journal of Medicine reported that among adults living in the United States, the levels of albumin in cerebrospinal fluid is normally distributed with a mean of 29.5 mg/dl, and a standard deviation of 9.25 mg/dl. We are going to select a sample of size 20.
a) How does the variability of our sample mean compare with the variability of albumin levels in the population?
Answer
The standard deviation of the sample mean is the standard error.
\(SE = \frac{\sigma}{\sqrt{n}} = \frac{9.25}{\sqrt{20}} = 2.07\)
This is LESS than the standard deviation of the albumin levels (9.25)
b) What is the probability that our sample mean will lie between 29 and 31 mg/dl?
Answer
When we are looking at the probability that our SAMPLE mean will lie betwen 2 values, we need to use the SE. Since we can assume the data are normally distributed we can still use the Z formula.
\[ P(29 < \bar{X} < 31) = P(\frac{29 - 29.5}{9.25/\sqrt{20}} < Z < \frac{31 - 29.5}{9.25/\sqrt{20}}) \\ = P(-0.242 < Z < 0.725) \\ = P(Z < 0.725) - P(Z < -0.242)\\ = 0.766 - 0.404 = 0.362 \]
c) What two values will contain the middle 50% of our sample means?
Answer
For the middle 50% of the sample means, we want to find \(z_{0.25}\) and \(z_{0.75}\) in the table. These 2 values are \(\pm0.675\). Since we are finding the middle 50% of the sample means rather than the individuals, we will use the standard error rather than the standard deviation. Therefore,
\(\mu \pm z *SE\) = \(29.5 \pm 0.675 * \frac{9.25}{\sqrt{20}}\) = \((28.104, 30.896)\)
d) Now assume we don’t know the standard deviation or mean of the population and our sample of 20 has a mean albumin level of 30.1 mg/dl and a standard deviation of 8.95 mg/dl. Construct a 95% confidence interval for the population mean and interpret it.
Answer
Since we no longer know the standard deviation of the population, we will use the T distribution:
\(95\% \space CI = \bar{x} \pm t_{0.025,
\space 19}{\frac{s}{\sqrt{n}}}\) = \(30.1 \pm 2.09*\frac{8.95}{\sqrt{20}}\) =
\((25.92, 34.28)\)
(note that if the course website table were used, if 19 degrees of freedom were on the table, 2.09 would be under \(\alpha = 0.05\))
We are 95% confident that the true mean level of albumnin in cerebrospinal fluid in adults living in the US is between 25.92 mg/dl and 34.28 mg/dl.
e) Why is the confidence interval constructed in (d) wider than the confidence interval that would be constructed if we used the normal distribution? What variables affect the width of the confidence interval?
Answer
Three values affect the width of the interval: the standard deviation, the sample size, and the critical value. In this case, since we are using the t-distribution, the tails are fatter and the critical value is larger, causing the interval to be wider.
6. In the following scenarios, identify what will happen to the power of a hypothesis test:
a) We increase the sample size.
Answer
Power increases
b) The standard deviation of the sample is larger than what we expected.
Answer
Power decreases
c) Our effect size moves from 5 units to 10 units.
Answer
Power increases