# A tibble: 3 × 4
Group mean variance sd
<fct> <dbl> <dbl> <dbl>
1 mAD 761. 182068. 427.
2 MCI 341. 165168. 406.
3 NCI 336. 189756. 436.
For some reason, these values get rounded to the whole number. You can get the full numbers to print if you use a package to print prettier tables, like kableExtra (only works for HTML output)
library(kableExtra)
Attaching package: 'kableExtra'
The following object is masked from 'package:dplyr':
group_rows
Call:
lm(formula = Abeta ~ Group, data = Amyloid, contrasts = list(Group = "contr.treatment"))
Residuals:
Min 1Q Median 3Q Max
-623.3 -322.3 -108.3 227.9 1269.0
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 761.3 102.4 7.432 8.19e-10 ***
GroupMCI -420.2 137.8 -3.050 0.00354 **
GroupNCI -425.0 141.0 -3.014 0.00392 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 422.3 on 54 degrees of freedom
Multiple R-squared: 0.1811, Adjusted R-squared: 0.1508
F-statistic: 5.971 on 2 and 54 DF, p-value: 0.004544
Intercept = mean for reference group
Mean for mAD = 761.3
GroupMCI is the difference between reference group and MCI
Mean for MCI = 761.3 - 420.2 = 341.1
GroupNCI is the difference between reference group and NCI
Mean for NCI = 761.3 - 425.0 = 336.3
Compare these values to the group means calculated above
t tests tell you whether that contrast (difference) is significant
Both contrasts are significant here
Mean Abeta for mAD is different from mean Abeta for MCI
Mean Abeta for mAD is different from mean Abeta for NCI
Report t tests as we’ve dont before
The degrees of freedom for each test is 1
5.1 Multiple comparisons
If we controlled for multiple comparisons what are our conclusions?
Do those differ from the conclusions where we didn’t control for multiple comparisons?
p.adjust(c(0.00354, 0.00392), "holm", 2)
[1] 0.00708 0.00708
p.adjust(c(0.00354, 0.00392), "BH", 2)
[1] 0.00392 0.00392
6 Conclusions
In plain language:
Are there differences in amyloid beta levels for the different impairment groups?
Which groups are different from one another?
Source Code
---title: "BTS 510 Lab 12"format: html: embed-resources: true self-contained-math: true html-math-method: katex number-sections: true toc: true code-tools: true code-block-bg: true code-block-border-left: "#31BAE9"---```{r}#| label: setuplibrary(tidyverse)set.seed(12345)theme_set(theme_classic(base_size =16))```## Learning objectives* *Describe* the logic of ANOVA* *Conduct* and *interpret* the overall ANOVA test* *Conduct* post hoc tests## Data* `Amyloid` dataset (N = 57) in **Stat2Data** package * `Group`: mAD=Alzheimer’s, MCI=mild impairment, NCI=no impairment * `Abeta`: Amount of Abeta from the posterior cingulate cortex (pmol/g tissue)```{r}library(Stat2Data)data(Amyloid)head(Amyloid)```## Descriptives and data visualization### Means, variances, and SDs by group```{r}Amyloid %>%group_by(Group) %>%summarize(mean =mean(Abeta),variance =var(Abeta),sd =sd(Abeta))```* For some reason, these values get rounded to the whole number. You can get the full numbers to print if you use a package to print prettier tables, like **kableExtra** (only works for HTML output) ```{r}library(kableExtra)Amyloid %>%group_by(Group) %>%summarize(mean =mean(Abeta),variance =var(Abeta),sd =sd(Abeta)) %>%kable()```### Equal variance assumption* Graphically```{r}ggplot(data = Amyloid,aes(x = Group, y = Abeta)) +geom_boxplot(alpha =0.5) +geom_jitter(width =0.2)```* Statistical test```{r}library(car)leveneTest(Abeta ~ Group, data = Amyloid)```* Variances are roughly equal -- good to go!## Analysis of variance* Use both `lm()` with `contr.sum` and `anova()` to get correct type III sums of squares```{r}m1 <-lm(Abeta ~ Group,contrasts=list(Group='contr.sum'),data = Amyloid)anova(m1)```* Overall, there is an effect of `Group` * We reject the $H_0$ that all the group means are equal, F(2, 54) = 5.97, p <.01. At least one group mean is different from the others.### Normality assumption* Graphically```{r}m1resid <-data.frame(residuals(m1))ggplot(data = m1resid, aes(sample = residuals.m1.)) +geom_qq() +geom_qq_line()m1resid <-data.frame(residuals(m1))ggplot(data = m1resid, aes(x = residuals.m1.)) +geom_histogram()```* Statistical test```{r}shapiro.test(resid(m1))```* Normality assumption is violated * We should consider non-parametric alternatives that don't assume normality * But for the purposes of this class, we're going to carry on anyway## Post hoc tests* I'm only going to use the reference / dummy / `treatment` contrasts here * Feel free to try out the others* Dummy coding compares each group to a specific reference group (e.g., control) * In R, the default reference group is the **first** group * Factors in R are ordered alphabetically (if character names) or numerically (if number names) * So this will do two comparisons (`mAD` is the reference group) * `mAD` vs `MCI` * `mAD` vs `NCI````{r}m1_trt <-lm(Abeta ~ Group,contrasts=list(Group='contr.treatment'),data = Amyloid)summary(m1_trt)```* Intercept = mean for reference group * Mean for `mAD` = 761.3* `GroupMCI` is the difference between reference group and `MCI` * Mean for `MCI` = 761.3 - 420.2 = 341.1* `GroupNCI` is the difference between reference group and `NCI` * Mean for `NCI` = 761.3 - 425.0 = 336.3* Compare these values to the group means calculated above* $t$ tests tell you whether that contrast (difference) is significant * Both contrasts are significant here * Mean `Abeta` for `mAD` is different from mean `Abeta` for `MCI` * Mean `Abeta` for `mAD` is different from mean `Abeta` for `NCI` * Report $t$ tests as we've dont before * The degrees of freedom for each test is 1### Multiple comparisons* If we controlled for multiple comparisons what are our conclusions? * Do those differ from the conclusions where we didn't control for multiple comparisons?```{r}p.adjust(c(0.00354, 0.00392), "holm", 2)p.adjust(c(0.00354, 0.00392), "BH", 2)```## ConclusionsIn plain language:1. Are there differences in amyloid beta levels for the different impairment groups?2. Which groups are different from one another?
