Homework 2: Solutions
Min Wang
October 13, 2026
1 Problem 1: Comparing two B12 measurement methods
Because both methods were used on each adult, the measurements are paired. Let \(D=\text{Method 1}-\text{Method 2}\). The null hypothesis is that the population mean (or, for the nonparametric paired tests, the center) of the within-person differences is zero. All tests below are two-sided.
method1 <- c(204, 238, 209, 277, 197, 226, 203, 131, 282, 76)
method2 <- c(199, 230, 198, 253, 180, 209, 213, 137, 250, 82)
d <- method1 - method2
b12 <- data.frame(Adult = 1:10, Method1 = method1, Method2 = method2,
Difference = d)
b12c(n = length(d), mean_difference = mean(d), sd_difference = sd(d))## n mean_difference sd_difference
## 10.00000 9.20000 13.798551.1 (a) Paired Student’s t-test and normality check
The paired t-test is a one-sample t-test on the differences. Its assumptions include independent adults and an approximately normal distribution of the within-adult differences. The Shapiro-Wilk test does not reject normality here (\(p\approx 0.747\)); the small sample means this check has limited power. The dot plot and normal Q-Q plot are also shown.
shapiro.test(d)##
## Shapiro-Wilk normality test
##
## data: d
## W = 0.95667, p-value = 0.7473par(mfrow = c(1, 1))
qqnorm(d, pch = 19, main = "Normal Q-Q plot of differences")
qqline(d, col = "steelblue", lwd = 2)
par(mfrow = c(1, 1))t.test(method1, method2, paired = TRUE, alternative = "two.sided",
conf.level = 0.99)##
## Paired t-test
##
## data: method1 and method2
## t = 2.1084, df = 9, p-value = 0.06424
## alternative hypothesis: true mean difference is not equal to 0
## 99 percent confidence interval:
## -4.980608 23.380608
## sample estimates:
## mean difference
## 9.2At the 1% level, the paired t-test does not reject equality of the method means: the two-sided p-value is about 0.064, which exceeds 0.01. The observed average difference is 9.2 units (Method 1 higher), but the evidence is not strong enough to conclude that the methods differ at this significance level.
1.2 (b) Sign test
Ignore zero differences (there are none). Under the null, positive and negative differences are equally likely. Seven of the ten differences are positive.
n_positive <- sum(d > 0)
n_nonzero <- sum(d != 0)
binom.test(n_positive, n_nonzero, p = 0.5, alternative = "two.sided")##
## Exact binomial test
##
## data: n_positive and n_nonzero
## number of successes = 7, number of trials = 10, p-value = 0.3438
## alternative hypothesis: true probability of success is not equal to 0.5
## 95 percent confidence interval:
## 0.3475471 0.9332605
## sample estimates:
## probability of success
## 0.7The exact two-sided p-value is 0.34375. Since this exceeds 0.01, the sign test does not reject the null hypothesis.
1.3 (c) Wilcoxon test
The Wilcoxon rank-sum test (Mann-Whitney test) is designed for two independent samples. These measurements are paired, so applying it to the two method columns discards the adult-level matching. It is included below to answer the requested calculation, but it is not the appropriate primary analysis.
wilcox.test(method1, method2, paired = FALSE, alternative = "two.sided",
exact = FALSE, correct = TRUE)##
## Wilcoxon rank sum test with continuity correction
##
## data: method1 and method2
## W = 54.5, p-value = 0.7623
## alternative hypothesis: true location shift is not equal to 0The rank-sum test also does not reject at 1% (its p-value is approximately 0.762). For paired data, the corresponding Wilcoxon procedure is the Wilcoxon signed-rank test, applied to the within-adult differences:
wilcox.test(d, mu = 0, alternative = "two.sided", exact = FALSE,
correct = TRUE)##
## Wilcoxon signed rank test with continuity correction
##
## data: d
## V = 45, p-value = 0.08273
## alternative hypothesis: true location is not equal to 0This paired signed-rank test does not reject at 1% either (approximate p-value about 0.083). The data contain tied absolute differences, so R uses the normal approximation rather than an exact signed-rank calculation.
1.4 (d) Comparison and practical choice
All four calculated tests fail to reject their respective null hypotheses at the specified 1% level. This common decision does not make the rank-sum test an appropriate analysis: the study design is paired. In practice, use the paired t-test when the within-person differences are reasonably normal; the observed differences show no clear departure from normality, and the paired t-test is the direct test of a mean difference. If a distribution-free paired analysis is preferred, use the sign test (few assumptions, but lower power) or the Wilcoxon signed-rank test (which uses both signs and ranks and assumes a symmetric distribution of differences). Do not use the rank-sum test for these matched observations.
2 Problem 2: Happiness and marital status
Let the row variable be happiness category and the column variable be marital status.
- \(H_0\): Happiness and marital status are independent in the population.
- \(H_a\): Happiness and marital status are associated in the population.
The Pearson chi-square test of independence uses significance level 0.05.
happiness <- matrix(c(600, 63, 112, 144,
720, 142, 355, 459,
93, 51, 119, 127),
nrow = 3, byrow = TRUE,
dimnames = list(
Happiness = c("Very Happy", "Pretty Happy", "Not Too Happy"),
MaritalStatus = c("Married", "Widowed", "Divorced/Separated",
"Never Married")))
happiness## MaritalStatus
## Happiness Married Widowed Divorced/Separated Never Married
## Very Happy 600 63 112 144
## Pretty Happy 720 142 355 459
## Not Too Happy 93 51 119 127chisq.test(happiness, correct = FALSE)##
## Pearson's Chi-squared test
##
## data: happiness
## X-squared = 224.12, df = 6, p-value < 2.2e-16The test statistic is \(\chi^2(6)=224.12\), with \(p<0.001\) (approximately \(1.38\times 10^{-45}\)). All expected cell counts exceed 5, so the usual large-sample condition is satisfied. At the 5% level, reject independence: there is strong evidence that happiness and marital status are associated. This is an association, not evidence that marital status causes happiness.
For context, the within-status proportions are:
round(prop.table(happiness, margin = 2), 3)## MaritalStatus
## Happiness Married Widowed Divorced/Separated Never Married
## Very Happy 0.425 0.246 0.191 0.197
## Pretty Happy 0.510 0.555 0.606 0.629
## Not Too Happy 0.066 0.199 0.203 0.174These proportions describe the pattern in this sample; they do not change the chi-square test’s conclusion.
3 Problem 3: Cholesterol by blood pressure status
Place heartbpchol.csv in the same folder as this R
Markdown file. The provided file contains the alive-patient records with
columns Cholesterol and BP_Status. The
supplied instructions say the normality assumption is reasonable. We fit
a one-way ANOVA and use Levene’s test to assess equal variances. The
ANOVA model’s \(R^2\) is the fraction
of observed cholesterol variation explained by blood pressure
status.
heart <- read.csv("C://Users/likeb/OneDrive/Desktop/STA6443/Homework/Homework 2/heartbpchol.csv", check.names = FALSE)
bp_col <- if ("BP_Status" %in% names(heart)) "BP_Status" else "BP Status"
needed <- c("Cholesterol", bp_col)
heart <- heart[complete.cases(heart[, needed]), needed]
names(heart)[names(heart) == bp_col] <- "BP_Status"
heart$BP_Status <- factor(heart$BP_Status, levels = c("Optimal", "Normal", "High"))
heart$BP_Status <- droplevels(heart$BP_Status)
fit_chol <- aov(Cholesterol ~ BP_Status, data = heart)
summary(fit_chol)## Df Sum Sq Mean Sq F value Pr(>F)
## BP_Status 2 25211 12605 6.671 0.00137 **
## Residuals 538 1016631 1890
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1model_summary <- summary(fit_chol)[[1]]
ss_model <- model_summary[1, "Sum Sq"]
ss_total <- sum(model_summary[, "Sum Sq"])
r2_chol <- ss_model / ss_total
cat(sprintf("R-squared = %.4f (%.1f%% of sample variation)\n",
r2_chol, 100 * r2_chol))## R-squared = 0.0242 (2.4% of sample variation)if (!requireNamespace("car", quietly = TRUE)) {
stop("Install package 'car' to run Levene's test: install.packages('car')")
}
car::leveneTest(Cholesterol ~ BP_Status, data = heart, center = median)aggregate(Cholesterol ~ BP_Status, data = heart,
FUN = function(z) c(n = length(z), mean = mean(z), sd = sd(z)))Interpret the output as follows: the ANOVA p-value tests whether all three population means are equal. If \(p<0.05\), conclude that at least one mean differs; if \(p\geq0.05\), the data do not provide sufficient evidence of a difference. Report the displayed \(R^2\) as the percentage of sample variability in cholesterol accounted for by BP Status. For Levene’s test, a small p-value (below 0.05) is evidence against equal population variances; a larger p-value means the equal-variance assumption is not contradicted by this check, but does not prove the variances are equal.
If the overall ANOVA is significant, Tukey’s HSD identifies which pairs differ while controlling the family-wise error rate:
TukeyHSD(fit_chol, "BP_Status", conf.level = 0.95)## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = Cholesterol ~ BP_Status, data = heart)
##
## $BP_Status
## diff lwr upr p adj
## Normal-Optimal 7.103198 -6.981749 21.18815 0.4624869
## High-Optimal 18.646679 4.456460 32.83690 0.0059898
## High-Normal 11.543481 2.153023 20.93394 0.0111929In the Tukey table, a pair is statistically different at 5% when its
adjusted p adj is below 0.05 (equivalently, its 95%
confidence interval does not contain zero). Compare group means above to
identify which status has the largest and smallest sample mean. For this
dataset, the ANOVA is significant, \(F(2,538)=6.67, p=0.0014\). BP Status
accounts for about 2.42% of the sample variation in cholesterol (\(R^2=0.0242\)), a statistically detectable
but modest share. Levene’s test is not significant (\(F=0.18, p=0.833\)); the equal-variance
assumption is reasonable based on this diagnostic. Tukey’s comparisons
show that the High mean (240.57) exceeds both the Normal mean (229.03;
adjusted \(p=0.011\)) and the Optimal
mean (221.93; adjusted \(p=0.006\)).
The Normal and Optimal means do not differ significantly (adjusted \(p=0.462\)). Thus, the High BP Status group
has the highest average cholesterol in this sample; this describes an
association and does not establish causation.
4 Problem 4: Blood tests by drink group
Place bupa.csv in the same folder. It must contain
numeric columns mcv and alkphos and a
categorical column drinkgroup. The assignment states that
normality and independence are reasonable. Each response is analyzed
with a separate one-way ANOVA; Levene’s test evaluates the
equal-variance assumption.
bupa <- read.csv("C://Users/likeb/OneDrive/Desktop/STA6443/Homework/Homework 2/bupa.csv", check.names = FALSE)
needed_bupa <- c("mcv", "alkphos", "drinkgroup")
if (!all(needed_bupa %in% names(bupa))) {
stop("bupa.csv must contain columns mcv, alkphos, and drinkgroup.")
}
bupa <- bupa[complete.cases(bupa[, needed_bupa]), needed_bupa]
bupa$drinkgroup <- factor(bupa$drinkgroup)
if (!requireNamespace("car", quietly = TRUE)) {
stop("Install package 'car' to run Levene's test: install.packages('car')")
}4.1 (a) MCV by drink group
fit_mcv <- aov(mcv ~ drinkgroup, data = bupa)
summary(fit_mcv)## Df Sum Sq Mean Sq F value Pr(>F)
## drinkgroup 4 733 183.29 10.26 7.43e-08 ***
## Residuals 340 6073 17.86
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tab_mcv <- summary(fit_mcv)[[1]]
r2_mcv <- tab_mcv[1, "Sum Sq"] / sum(tab_mcv[, "Sum Sq"])
cat(sprintf("R-squared = %.4f (%.1f%% of sample variation)\n",
r2_mcv, 100 * r2_mcv))## R-squared = 0.1077 (10.8% of sample variation)car::leveneTest(mcv ~ drinkgroup, data = bupa, center = median)aggregate(mcv ~ drinkgroup, data = bupa,
FUN = function(z) c(n = length(z), mean = mean(z), sd = sd(z)))TukeyHSD(fit_mcv, "drinkgroup", conf.level = 0.95)## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = mcv ~ drinkgroup, data = bupa)
##
## $drinkgroup
## diff lwr upr p adj
## 2-1 1.241452991 -0.690316778 3.173223 0.3973587
## 3-1 0.938131313 -0.697363908 2.573627 0.5157202
## 4-1 3.744610282 1.968846495 5.520374 0.0000002
## 5-1 3.746031746 0.999127120 6.492936 0.0020039
## 3-2 -0.303321678 -2.330666517 1.724023 0.9940309
## 4-2 2.503157290 0.361050967 4.645264 0.0127884
## 5-2 2.504578755 -0.492214115 5.501372 0.1499436
## 4-3 2.806478969 0.927189295 4.685769 0.0005031
## 5-3 2.807900433 -0.007037875 5.622839 0.0509321
## 5-4 0.001421464 -2.897262735 2.900106 1.0000000The ANOVA p-value determines whether there is evidence that at least one population mean MCV differs among drink groups. The \(R^2\) reports the fraction of sample MCV variation explained by group. Levene’s p-value below 0.05 suggests unequal variances; otherwise equal variance is not rejected by this diagnostic. In the Tukey output, adjusted p-values below 0.05 identify pairs with different means. Use the group means to state which groups have the highest and lowest observed MCV.
4.2 (b) Alkaline phosphatase by drink group
fit_alkphos <- aov(alkphos ~ drinkgroup, data = bupa)
summary(fit_alkphos)## Df Sum Sq Mean Sq F value Pr(>F)
## drinkgroup 4 4946 1236.4 3.792 0.00495 **
## Residuals 340 110858 326.1
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tab_alkphos <- summary(fit_alkphos)[[1]]
r2_alkphos <- tab_alkphos[1, "Sum Sq"] / sum(tab_alkphos[, "Sum Sq"])
cat(sprintf("R-squared = %.4f (%.1f%% of sample variation)\n",
r2_alkphos, 100 * r2_alkphos))## R-squared = 0.0427 (4.3% of sample variation)car::leveneTest(alkphos ~ drinkgroup, data = bupa, center = median)aggregate(alkphos ~ drinkgroup, data = bupa,
FUN = function(z) c(n = length(z), mean = mean(z), sd = sd(z)))TukeyHSD(fit_alkphos, "drinkgroup", conf.level = 0.95)## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = alkphos ~ drinkgroup, data = bupa)
##
## $drinkgroup
## diff lwr upr p adj
## 2-1 -2.645299 -10.8987258 5.608127 0.9045115
## 3-1 -4.056138 -11.0437411 2.931464 0.5035761
## 4-1 -1.148743 -8.7356393 6.438152 0.9937509
## 5-1 12.572650 0.8365845 24.308715 0.0288892
## 3-2 -1.410839 -10.0726073 7.250929 0.9917361
## 4-2 1.496556 -7.6555274 10.648639 0.9916117
## 5-2 15.217949 2.4142438 28.021654 0.0107154
## 4-3 2.907395 -5.1218122 10.936602 0.8583931
## 5-3 16.628788 4.6020510 28.655525 0.0016498
## 5-4 13.721393 1.3368544 26.105932 0.0214375Interpret the ANOVA, \(R^2\), Levene, and Tukey results in the same way as for MCV. A significant overall ANOVA indicates that at least one group mean differs, while Tukey’s HSD shows which pairs differ. The group means identify the highest and lowest observed alkaline phosphatase levels.
4.3 (c) Comparing the two outcomes
Compare the Tukey tables and the group means for the two outcomes.
Note which pairwise group differences are significant for MCV, which are
significant for alkaline phosphatase, and whether the same drink groups
tend to have higher or lower means for both tests. The two outcomes need
not show the same pattern: MCV and alkaline phosphatase measure
different blood characteristics. Avoid claiming that one test is more
sensitive based only on p-values; also consider the estimated
differences and their confidence intervals. The bupa.csv
file was not included with the homework PDF, so these numerical
conclusions are filled in by running the code with the course data
file.