\(\chi^2\)-Test

Wir stellen uns am Anfang eine ganz wichtige Frage! Ist die Ernährungsweise unabhängig vom Geschlecht? Ja, diese Frage ist wichtig. Das wichtige an ihr für uns ist vor allem, dass sie sich wunderbar für einen \(\chi^2\)-Test eignet.

result <- chisq.test(data$geschlecht, data$vegetarier)
result

## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  data$geschlecht and data$vegetarier
## X-squared = 3.9077, df = 1, p-value = 0.04807

result$observed

##                data$vegetarier
## data$geschlecht ja nein
##        männlich  2   21
##        weiblich 27   57

Wie zu erwarten: Männer wollen immer noch ihr Steak. Der Effekt ist gar nicht klein und auf einem \(\alpha\)-Niveau von 5% haben wir ein signifikantes Ergebnis. Aufpassen müssen wir, wenn wir das Ergebnis per Hand nachrechnen wollen, denn man sieht im Output, dass hier noch eine Korrektur verwendet wird, die zu einem genaueren \(p\)-Wert führt. Das Problem ist, dass sich dadurch der \(\chi^2\)-Wert verändert und der Effekt würde falsch berechnet werden:

sqrt(3.9077 / 107) # für Effektgröße w und bei 2x2-Tabelle ist das auch Phi

## [1] 0.1911035

Tatsächlich ist der Effekt:

cor(as.numeric(as.factor(data$geschlecht)), as.numeric(as.factor(data$vegetarier)))

## [1] -0.2166947

Wobei hier die spezifische Kodierung zu einem negativen Effekt führt, was für uns allerdings keine Rolle spielt, denn wir sehen aus der Tabelle, dass Frauen eher Vegetarier sind als Männer. Die komischen as.factor und as.numeric sind hier notwendig, denn das Geschlecht und die Ernährungsweise sind als Characters gespeichert.

Wie kommen wir über den \(\chi^2\)-Wert auf den richtigen Effekt? Einfach die Korrektur abschalten:

chisq.test(data$geschlecht, data$vegetarier, correct = FALSE)

## 
##  Pearson's Chi-squared test
## 
## data:  data$geschlecht and data$vegetarier
## X-squared = 5.0244, df = 1, p-value = 0.02499

Und nun passt die Formel:

sqrt(5.0244/107)

## [1] 0.2166957

Die Detektivarbeit war erfolgreich und wir sehen auch, dass die Korrelation äquivalent zur Effektgröße \(\Phi\) und \(w\) ist – zumindest für 2×2-Tabellen. Jetzt bist Du dran. In der Übung zur Methodenlehre haben wir aus einem Artikel folgende Tabelle extrahiert:

a <- c(19, 27, 9, 13, 8, 27)
b <- c(25, 26, 24, 15, 9, 4)
gutachter <- matrix(a+b, nrow=2, byrow=T)
colnames(gutachter) <- c("Plant Relocation", "Tax evasion", "Independence")
gutachter <- as.data.frame(gutachter, row.names = c("Recognize", "Do not recognize"))
gutachter

Es ging um die Frage ob Gutachter ethische Probleme in Finanz-Anträgen identifizieren können. Es gab drei Szenarien (in den Spalten) und als AV wurde erfasst, ob das Problem erkannt wurde oder nicht. Man kann schon nach einem kurzen Blick erkennen, dass bei Steuerhinterziehung die Erkennungsrate am höchsten ist. Die Frage ist generell, ob es Unterschiede in der Erkennungsrate zwischen den Szenarien gibt. Als \(\chi^2\)-Wert kam im Paper 8.6 heraus. Prüfe diesen nach:

chisq.test(gutachter)

\(U\)-Test

Wir haben an Tag 4 einen \(t\)-Test gerechnet um zu prüfen, ob sich die Körpergröße zwischen Männern und Frauen unterscheidet. Die Voraussetzungen für diesen Test sind vielleicht nicht erfüllt. Eine Alternative ist hier der \(U\)-Test, der in R jedoch unter dem Wilcoxon-Test läuft (es gibt viele Namen für diesen Test, auch noch Mann-Whitney-Wilcoxon, o. Wilcoxon-Mann-Whitney o. eben Wilcoxon rank-sum test):

m <- data[data$geschlecht == "männlich", "groesse"]
w <- data[data$geschlecht == "weiblich", "groesse"]
wilcox.test(m, w, paired = FALSE) # paired = FALSE ist der default, kann man also auch weglassen

## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  m and w
## W = 1793.5, p-value = 3.403e-10
## alternative hypothesis: true location shift is not equal to 0

Im Output sieht man mal wieder, dass eine Korrektur angewandt wird. Wenn Du also was per Hand rechnest, vergiss nicht die Korrektur auszustellen. Übrigens ist \(W\) hier die Teststatistik, die wir eigentlich als \(U\) kennen.

Ist gar nicht schwer, oder? Das kannst Du auch. In der Übung haben wir uns ein Paper angeschaut, bei dem es um die Frage ging ob Frauen und Männer in unterschiedlichen Umwelten (friedlich oder kompetitiv) besser oder schlechter performen. Im Folgenden findest Du die Anzahl der gelösten Aufgaben:

solved_female <- c(24, 20, 18, 17, 15, 15, 14, 14, 14, 13, 13, 11, 11, 10, 9, 
                   9, 8, 8, 8, 8, 8, 8, 8, 7, 7, 7, 7, 7, 4, 4)
solved_male <- c(22, 21, 21, 21, 19, 19, 17, 17, 17, 17, 17, 18, 16, 16, 16,
                 15, 15, 15, 15, 13, 13, 13, 12, 11, 11, 11, 10, 10, 8, 5)

Rechne den \(U\)-Test für den Vergleich zwischen Männern und Frauen. Stell hier die Korrektur aus, sodass der gleiche Wert herauskommt wie bei der Berechnung per Hand in der Übung.

wilcox.test(solved_female, solved_male, correct = FALSE)

Wilcoxon-Test (signed-rank)

Wir haben an Tag 4 auch einen \(t\)-Test gerechnet um zu prüfen, ob sich die Deutsch- und Mathepunkte unterscheiden. Die Voraussetzungen für diesen Test sind jedoch sehr wahrscheinlich nicht erfüllt. Eine Alternative ist hier der Wilcoxon-Test (signed-rank):

wilcox.test(data$mathe, data$deutsch, paired = TRUE)

## 
##  Wilcoxon signed rank test with continuity correction
## 
## data:  data$mathe and data$deutsch
## V = 1773, p-value = 0.202
## alternative hypothesis: true location shift is not equal to 0

Die Test-Statistik ist hier \(V\), wir kennen sie unter \(T\).

Die entscheidende Frage ist, wie weit der \(p\)-Wert dieses Tests vom \(p\)-Wert des \(t\)-Tests entfernt ist:

t.test(data$mathe, data$deutsch, paired = TRUE)

## 
##  Paired t-test
## 
## data:  data$mathe and data$deutsch
## t = -1.623, df = 106, p-value = 0.1076
## alternative hypothesis: true mean difference is not equal to 0
## 95 percent confidence interval:
##  -0.9965795  0.0993832
## sample estimates:
## mean difference 
##      -0.4485981

Tatsächlich macht das hier einen großen Unterschied (0.202 versus 0.1076)! Das ist jedoch nicht immer so, man sollte sich aber auf jeden Fall Gedanken machen auf welchem Skalenniveau man misst und ob die Voraussetzungen für einen Test näherungsweise gegeben sind. Ist dies nicht der Fall, so ist man mit den nonparametrischen Tests auf der sicheren Seite. Wie schon erwähnt, die eigentlich bessere Alternative ist Bootstrapping und Randomisierungstests, aber die behandeln wir erst im Master.

Interessant wäre hier noch zu schauen, wie groß der Unterschied zu einem unabhängigen Wilcoxon-Test ist:

wilcox.test(data$mathe, data$deutsch, paired = FALSE)

## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data$mathe and data$deutsch
## W = 5396, p-value = 0.4648
## alternative hypothesis: true location shift is not equal to 0

Macht einen noch größeren Unterschied (0.202 versus 0.4648)! Das sollte klar sein, denn bei abhängigen Stichproben sind die Effekte immer größer und man sollte das hier nicht verschenken. Im umgekehrten Fall sollte man auch nie unabhängige Stichproben als abhängige behandeln; das wäre grob fahrlässig!

In der Übung hatten wir auch einen Wilcoxon-Test gerechnet. Es ging um eine neue Therapieform bei Anorexie-Patienten. Gemessen wurde das Körpergewicht in kg vor und nach der Therapie. Die Daten waren:

gewicht_t1 <- c(58, 55, 60, 52, 53, 53, 55, 49, 50, 55)
gewicht_t2 <- c(61, 60, 64, 56, 59, 60, 59, 47, 52, 56)

Prüfe nach ob Du auf das selbe Ergebnis kommst wie in der Übung (wir haben einseitig getestet, auf eine Erhöhung des Gewichts):

wilcox.test(gewicht_t2, gewicht_t1, paired = TRUE, alternative="greater", correct = FALSE)

Sehr gut! Jetzt haben wir aber wirklich alle Signifikanztests der Vorlesung behandelt. Ich hoffe ich habe nichts vergessen.

Das war harte Arbeit und jetzt kommt die Belohnung. Poweranalyse :). Klingt wohl eher nach noch mehr Arbeit, aber bevor Du zu nicht-freier, unmoralischer Software wie G*Power greifst, investiere lieber noch 30 Minuten in R!

Aufgabe \(\chi^2\)-Test

Die folgenden Daten sind aus einer Übungsaufgabe der Methodenlehre zum Simpsons Paradox.

a <- c(820, 80, 680, 20)
b <- c(20, 80, 100, 200)
simpsonsparadox <- matrix(a+b, nrow=2, byrow=T)
colnames(simpsonsparadox) <- c("Accept", "Reject")
simpsonsparadox <- as.data.frame(simpsonsparadox, row.names = c("Male", "Female"))
simpsonsparadox

Dargestellt sind hier die Häufigkeiten der Annahme und Ablehnung von Studienplatzbewerbungen aufgeteilt nach Geschlecht. Wir haben in der Übung den \(p\)-Wert aus dem Artikel von Kievit hierzu nachgerechnet. Das haben wir per Hand gemacht, also ohne Yates Korrektur und kamen auf einen Wert von 11.696. Welcher Wert kommt in R (mit Korrektur) heraus?

chisq.test(simpsonsparadox, correct = T)

Das ist auch der Wert, der im Paper angegeben ist.

Aufgabe Partialkorrelation

Schreibe eine Funktion mit der Du aus Korrelationen Partialkorrelationen berechnen kannst. Benutze die Funktion anschließend für das Beispiel aus der Vorlesung:

vl <- matrix(c(1, .43, -.51, -.54, .43, 1, .41, -.12, -.51, .41, 1, .23, -.54, -.12, .23, 1), ncol = 4)
row.names(vl) <- c("Gewichtsverlust", "Sport", "Kalorienaufnahme", "Alter")
colnames(vl) <- c("Gewichtsverlust", "Sport", "Kalorienaufnahme", "Alter")
vl 

##                  Gewichtsverlust Sport Kalorienaufnahme Alter
## Gewichtsverlust             1.00  0.43            -0.51 -0.54
## Sport                       0.43  1.00             0.41 -0.12
## Kalorienaufnahme           -0.51  0.41             1.00  0.23
## Alter                      -0.54 -0.12             0.23  1.00

Der Zusammenhang zwischen Gewichtsverlust und Sport, ohne den Einfluss von Kalorienaufnahme soll berechnet werden. Achtung: Das ist für Anfänger eine recht anspruchsvolle Aufgabe, für die Du etwas Zeit einplanen solltest.

parcor <- function(xy, xz, yz){
  (xy - xz * yz) / (sqrt(1 - xz^2) * sqrt(1 - yz^2))
}
parcor(vl[1, 2], vl[1, 3], vl[2, 3])

Aufgabe Power

Wie groß ist die Power bei einer Regression mit zwei Prädiktoren, einer Stichprobengröße von 88, einem \(\alpha\) von 5% und einem \(R^2\) von 0.13? Benutze die Funktion wp.regression und finde heraus, wie Du \(R^2\) in \(f^2\) umrechnen kannst.

f2 = 0.13 / (1 - 0.13)
wp.regression(n = 88, p1 = 2, f2 = f2, alpha = .05)

Gut gemacht. Du hast diesen Tag erfolgreich überstanden. Nur noch zwei Tage, dann bist Du ein R-Profi!