library(Rcmdr)
data("ChickWeight") # load data set
Nachschlagen, was die Spalten in den Daten bedeuten
#?ChickWeight
Unten ist bereits ein Code, der eine Abbildung erzeugt, die den Zusammenhang zwischen Zeit und Gewicht für jedes Küken abbildet.
Kopieren Sie den Code, um eine zweite Abbildung zu erstellen. Ändern Sie den Code dabei so ab, dass der Zusammehang zwischen Zeit und Gewicht für die verschiedenen Diäten sichtbar wird. Verändern Sie den Code auch so, dass die Grafik nur eine Zeile (und 4 Spalten, also eine pro Diät) hat. Färben Sie die Grafik dunkelgrün.
Erstellen Sie danach einen Strip-Chart, in dem lediglich die Gewichte der Küken für die verschiedenen Arten von Diät abgebildet werden. Alle vier Strips sollen verschiedene Farben haben. Für die Einfärbung sollen Sie den RColorBrewer verwenden. Wie das geht, können Sie hier lernen. Verwenden Sie für die Grafik die Farbpalette ‘RdBu’.
Erstellen Sie einen Plot für Mittelwerte und Konfidenzintervalle (vernachlässigen Sie hier den Faktor “Zeit”).
Weiter unten im Code ist bereits ein Zusat-Plot eingefügt, der mit dem Paket ggplot2 erstellt wurde. Das ist ein Paket, mit dem man (theoretisch) professionelle “Hochglanzabbildungen” erstellen kann. Wenn Sie mehr über das erstellen von Grafiken mit ggplot lernen möchten, schauen Sie sich dieses frei verfügbare Online-Buch an.
Erstellen Sie ein Subset des Datensatzes, das für die Zeitvariable nur den ersten und den letzten Tag der Messungen enthält.
Konvertieren Sie die Zeitvariable in einen Faktor. Vergeben Sie Labels für die Faktorstufen (“first day” und “last day”).
Machen Sie dann wieder eine Mittelwertsabbildung, die die Körpergewichte zum ersten und letzten Tag zeigt und dabei auch die vier verschiedenen Diäten berücksichtigt. Ändern Sie den Code für die Grafik so ab, dass Zeit auf der X-Achse ist und “Diet” die verschiedenen Linienfarben sind.
Erstellen Sie nun zwei Subsets: Eins mit “first day” und eins mit “last day”.
Prüfen Sie mittels einfaktorieller ANOVA die Hypothese, dass die Küken in den verschiedenen Diätgruppen an Tag 1 ähnliches Gewicht haben. Levene-Test nicht vergessen.
Prüfen Sie mittels einfaktorieller ANOVA die Hypothese, dass die Küken in den verschiedenen Diätgruppen am letzten Tag unterschiedliche Gewichte haben. Levene-Test nicht vergessen.
Zusammenhang von Zeit und Gewicht für jedes Küken. Der Code für die Abbildung ist in der R-Beschreibung für den Datensatz zu finden.
require(graphics)
coplot(weight ~ Time | Chick, data = ChickWeight,
type = "b", show.given = FALSE)
Abändern:
require(graphics)
coplot(weight ~ Time | Diet, data = ChickWeight,
type = "b", show.given = FALSE, columns = 4, col = "darkgreen")
Strip-Chart:
library("RColorBrewer")
stripchart(weight ~ Diet, vertical=TRUE, method="jitter", jitter = 0.25, xlab="Diet",
ylab="weight", col=brewer.pal(n = 4, name = "RdBu"), data=ChickWeight)
Mittelwerte nach Diäten:
with(ChickWeight, plotMeans(weight, Diet, error.bars="conf.int", level=0.95,
connect=TRUE))
### Zusatzgrafik mit ggplot2
library(ggplot2) # das Paket ggplot2 müssen Sie sich ggf. installieren
myTheme <- theme(plot.title = element_text(face="bold", size = 22),
axis.title.x = element_text(size = 20),
axis.title.y = element_text(size = 20),
axis.text.x = element_text(size = 16, angle = 0),
axis.text.y = element_text(size = 16, angle = 0),
legend.text = element_text(size = 18),
legend.title = element_text(face = "bold", size = 18),
strip.text.x = element_text(size = 18),
#panel.grid.major = element_blank(),
#panel.grid.minor = element_blank(),
panel.background = element_blank(),
panel.grid.major = element_line(color = "#bdbdbd", size = 0.5,linetype = 2),
axis.line.x = element_line(colour = "black"),
axis.line.y = element_line(colour = "black"),
axis.text = element_text(colour ="black"),
axis.ticks = element_line(colour ="black"))
data_graph <- ChickWeight
data_graph$Diet <- factor(data_graph$Diet, levels = c("1", "2", "3", "4"),
labels = c("Diet 1", "Diet 2", "Diet 3", "Diet 4"))
pd <- position_dodge(width = 0.3)
# new labes for the facets
g <- ggplot(data_graph, aes(x=Time, y=weight, group = Chick)) +
guides(fill=FALSE)+
facet_grid( ~ Diet)+
ggtitle("Trajectory of chick weight over time given \ndifferent diets") +
scale_y_continuous(limits = c(-5, 400), breaks=seq(0, 400, 50), expand = c(0,0)) +
geom_line(position = pd, color = "black", size = 1, alpha=0.2) +
geom_jitter(aes(color = Diet), alpha = 0.3, width = 0.1) +
stat_summary(aes(y = weight,group=1, fill = Diet), fun.data = mean_cl_boot, geom = "ribbon", width = 0,
size = 1, alpha = 0.5) +
stat_summary(aes(y = weight,group=1, col = Diet), fun.y=mean,
geom="line",group=1, size = 1.5, linetype = "solid", alpha = 1)+
#stat_summary(aes(y = weight,group=1, fill = Diet), fun.y=mean, geom="point",
# color = "black", shape = 22, size = 3, group=1, alpha = 0.8)+
#stat_summary(aes(y = weight,group=1), fun.y=median, geom="point",
# color = "black", shape = 3, size = 4, group=1, alpha = 1,
# position = position_dodge(width = 0.5))+
labs(x = "Time (in Days)", y = "Chick weight (g)") +
scale_color_manual(name = "Diet",values=brewer.pal(n = 4, name = "RdBu"))+
scale_fill_manual(name = "Diet",values=brewer.pal(n = 4, name = "RdBu"))+
theme(legend.position = "none")+
myTheme
## Warning: Ignoring unknown parameters: width
## Warning: `fun.y` is deprecated. Use `fun` instead.
#theme_bw()
g
#ggsave("results_means.svg",width=9,height=7) # mit diesen Codes können Grafiken gespeichert werden
#ggsave("results_means.pdf",width=9,height=9) # speichern als pdf
Daten subsetten und in Faktor konvertieren:
dat_1st_last <- subset(ChickWeight, Time == 0 | Time == 21) # Subset erstellen
dat_1st_last$Time <- factor(dat_1st_last$Time, levels = c(0,21), labels = c("first day", "last day")) # faktorisieren
Abbildung:
with(dat_1st_last, plotMeans(weight, Time, Diet, error.bars="conf.int",
level=0.95, connect=TRUE, legend.pos="farright"))
Subsets für first und last day
first_day <- subset(dat_1st_last, Time == "first day")
last_day <- subset(dat_1st_last, Time == "last day")
Mit einfakt. ANOVA testen, ob sich die Gewichte zu Tag 1 zwischen den Diäten unterscheiden. Zuvor aber Leven-Test, ob Varianzen zu Tag 1 gleich sind.
library(mvtnorm, pos=19)
library(survival, pos=19)
library(MASS, pos=19)
library(TH.data, pos=19)
##
## Attaching package: 'TH.data'
## The following object is masked from 'package:MASS':
##
## geyser
library(multcomp, pos=19)
library(abind, pos=24)
Tapply(weight ~ Diet, var, na.action=na.omit, data=first_day) # variances by group
## 1 2 3 4
## 0.9894737 2.2333333 1.0666667 1.1111111
leveneTest(weight ~ Diet, data=first_day, center="mean") # levene test
## Levene's Test for Homogeneity of Variance (center = "mean")
## Df F value Pr(>F)
## group 3 1.9615 0.133
## 46
Test nicht signifikant, also H0, dass Varianzen zu Tag 1 gleich sind, beibehalten. Jetzt kan “normale” ANOVA ohne Welch-Korrektur (für ungleiche Varianzen) durchgeführt werden.
AnovaModel.1 <- aov(weight ~ Diet, data=first_day) # erstellt ANOVA Modell
summary(AnovaModel.1) # Gibt das Ergebnis zurück
## Df Sum Sq Mean Sq F value Pr(>F)
## Diet 3 4.32 1.440 1.132 0.346
## Residuals 46 58.50 1.272
with(first_day, numSummary(weight, groups=Diet, statistics=c("mean", "sd"))) # Ab hier alles Plot für Mittel-Diffs
## mean sd data:n
## 1 41.4 0.9947229 20
## 2 40.7 1.4944341 10
## 3 40.8 1.0327956 10
## 4 41.0 1.0540926 10
local({
.Pairs <- glht(AnovaModel.1, linfct = mcp(Diet = "Tukey"))
print(summary(.Pairs)) # pairwise tests
print(confint(.Pairs, level=0.95)) # confidence intervals
print(cld(.Pairs, level=0.05)) # compact letter display
old.oma <- par(oma=c(0, 5, 0, 0))
plot(confint(.Pairs))
par(old.oma)
})
##
## Simultaneous Tests for General Linear Hypotheses
##
## Multiple Comparisons of Means: Tukey Contrasts
##
##
## Fit: aov(formula = weight ~ Diet, data = first_day)
##
## Linear Hypotheses:
## Estimate Std. Error t value Pr(>|t|)
## 2 - 1 == 0 -0.7000 0.4368 -1.603 0.385
## 3 - 1 == 0 -0.6000 0.4368 -1.374 0.519
## 4 - 1 == 0 -0.4000 0.4368 -0.916 0.795
## 3 - 2 == 0 0.1000 0.5043 0.198 0.997
## 4 - 2 == 0 0.3000 0.5043 0.595 0.932
## 4 - 3 == 0 0.2000 0.5043 0.397 0.978
## (Adjusted p values reported -- single-step method)
##
##
## Simultaneous Confidence Intervals
##
## Multiple Comparisons of Means: Tukey Contrasts
##
##
## Fit: aov(formula = weight ~ Diet, data = first_day)
##
## Quantile = 2.6601
## 95% family-wise confidence level
##
##
## Linear Hypotheses:
## Estimate lwr upr
## 2 - 1 == 0 -0.7000 -1.8618 0.4618
## 3 - 1 == 0 -0.6000 -1.7618 0.5618
## 4 - 1 == 0 -0.4000 -1.5618 0.7618
## 3 - 2 == 0 0.1000 -1.2416 1.4416
## 4 - 2 == 0 0.3000 -1.0416 1.6416
## 4 - 3 == 0 0.2000 -1.1416 1.5416
##
## 1 2 3 4
## "a" "a" "a" "a"
Die ANOVA mit dem between-subject-Faktor “Diet” ist nicht signifikant, \(F(3,46)~= 1.132\), \(p~= 346\). Die Hypothese, dass sich die Gewichte der Küken zu Tag 1 zwischen den vier verschiedenen Diäten nicht unterscheiden, wird beibehalten.
Jetzt Testen, ob die Gewichte am letzten Tag unterschiedlich sind (aber erste Levene-Test):
Tapply(weight ~ Diet, var, na.action=na.omit, data=last_day) # variances by group
## 1 2 3 4
## 3445.933 6105.567 5129.789 1879.028
leveneTest(weight ~ Diet, data=last_day, center="mean")
## Levene's Test for Homogeneity of Variance (center = "mean")
## Df F value Pr(>F)
## group 3 1.2412 0.3072
## 41
Varianzhomogenität kann angenommen werden.
Jetzt ANOVA:
AnovaModel.2 <- aov(weight ~ Diet, data=last_day)
summary(AnovaModel.2)
## Df Sum Sq Mean Sq F value Pr(>F)
## Diet 3 57164 19055 4.655 0.00686 **
## Residuals 41 167839 4094
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
with(last_day, numSummary(weight, groups=Diet, statistics=c("mean", "sd")))
## mean sd data:n
## 1 177.7500 58.70207 16
## 2 214.7000 78.13813 10
## 3 270.3000 71.62254 10
## 4 238.5556 43.34775 9
local({
.Pairs <- glht(AnovaModel.2, linfct = mcp(Diet = "Tukey"))
print(summary(.Pairs)) # pairwise tests
print(confint(.Pairs, level=0.95)) # confidence intervals
print(cld(.Pairs, level=0.05)) # compact letter display
old.oma <- par(oma=c(0, 5, 0, 0))
plot(confint(.Pairs))
par(old.oma)
})
##
## Simultaneous Tests for General Linear Hypotheses
##
## Multiple Comparisons of Means: Tukey Contrasts
##
##
## Fit: aov(formula = weight ~ Diet, data = last_day)
##
## Linear Hypotheses:
## Estimate Std. Error t value Pr(>|t|)
## 2 - 1 == 0 36.95 25.79 1.433 0.4854
## 3 - 1 == 0 92.55 25.79 3.588 0.0048 **
## 4 - 1 == 0 60.81 26.66 2.281 0.1186
## 3 - 2 == 0 55.60 28.61 1.943 0.2252
## 4 - 2 == 0 23.86 29.40 0.811 0.8480
## 4 - 3 == 0 -31.74 29.40 -1.080 0.7025
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Adjusted p values reported -- single-step method)
##
##
## Simultaneous Confidence Intervals
##
## Multiple Comparisons of Means: Tukey Contrasts
##
##
## Fit: aov(formula = weight ~ Diet, data = last_day)
##
## Quantile = 2.6762
## 95% family-wise confidence level
##
##
## Linear Hypotheses:
## Estimate lwr upr
## 2 - 1 == 0 36.9500 -32.0748 105.9748
## 3 - 1 == 0 92.5500 23.5252 161.5748
## 4 - 1 == 0 60.8056 -10.5401 132.1512
## 3 - 2 == 0 55.6000 -20.9761 132.1761
## 4 - 2 == 0 23.8556 -54.8190 102.5301
## 4 - 3 == 0 -31.7444 -110.4190 46.9301
##
## 1 2 3 4
## "a" "ab" "b" "ab"
Effektgröße ausrechnen:
library(lsr)
etaSquared(AnovaModel.2)
## eta.sq eta.sq.part
## Diet 0.2540591 0.2540591
Die ANOVA ist signifikant, \(F(3,41)~= 4.655\), \(p~= .007\), \(\eta^2~= 0.25\). Die Hypothese, dass sich die Gewichte am letzten Tag zwischen den vier verschiedenen Diäten unterscheiden, wird angenommen. Die post-hoc Tests zeigen, dass der Unterschied durch den Unterschied zw. Diät-1 und Diät-3 entsteht.