Kapitel 12 Deskriptive Statistik
In diesem Abschnitt widmen wir uns univariaten und bivariaten deskriptiven Statistiken, vom Mittelwert und Standardabweichung über \(\chi^2\), Somers’ D und Korrelationskoeffizienten.
Wir halten uns (so grob) an die Struktur aus QM1 und QM2, das heißt: Wir kümmern uns beispielsweise zuerst um \(\chi^2\) als deskriptive Statistik um dann auf Cramer’s V zu sprechen zu kommen, und widmen uns erst im folgenden Kapitel dem \(\chi^2\)-Test als inferenzstatistisches Werkzeug.
12.1 Univariate Statistiken
Univariat ist fancy sprech für “nur eine Variable”, das heißt es geht erstmal noch nicht um Zusammenhänge irgendeiner Art, sondern primär darum wie ein Merkmal verteilt ist.
Die dazugehörigen Kennwerte und Formeln sind sozusagen das kleine 1x1 der Statistik, dementsprechend sind sie auch in R schnell und einfach erreichbar.
Als Beispieldatensatz nehmen wir wieder den Game of Thrones Deaths-Datensatz und den qmsurvey
-Datensatz, beide aus dem Kapitel Datenimport:
library(readr)
qmsurvey <- readRDS("data/qm_survey_ss2017.rds")
gotdeaths <- read_csv("data/got_deaths.csv", col_types = cols())
12.1.1 Mittelwert, Median, Modus
Mittelwert und Median sind einfach:
#> [1] 3.539216
#> [1] 4
Die Funktionen mean
und median
sind da wenig überraschend.
Aber was machen wir mit dem Modus? In R an sich gibt es keine modus
-Funktion, deswegen haben wir an dieser Stelle zwei Optionen:
- Häufigkeitstabelle erstellen und ablesen.
- Ein package benutzen, dass eine Funktion für den Modus hat.
Eine Häufigkeitstabelle bekommen wir sehr einfach mit table
:
#>
#> 1 2 3 4 5
#> 2 13 25 52 10
Das ist zwar optisch wenig beeindruckend, aber wir können immerhin einfach erkennen, dass der Wert 4
die größte absolute Höufigkeit hat. Passt.
Alternativ können wir die tadaatoolbox
benutzen:
#> [1] "4"
Vollkommen unerwartet kommen wir zum gleichen Ergebnis.
Who’da thunk.
Für die EnthusiastInnen:
Die modus
-Funktion macht übrigens auch nichts anderes als eine Häufigkeitstabelle via table
, aus der dann das Maximum extrahiert wird.
Praktisch nur names(table(x)[table(x) == max(table(x))])
mit ein bisschen Glitzer und so.
Wenn wir das jetzt noch in schön und übersichtlich haben wollen, können wir auf sjmisc
zurückgreifen:
val | frq | raw.prc | valid.prc | cum.prc | |
---|---|---|---|---|---|
1 | 2 | 1.96 | 1.96 | 1.96 | |
2 | 13 | 12.75 | 12.75 | 14.71 | |
3 | 25 | 24.51 | 24.51 | 39.22 | |
4 | 52 | 50.98 | 50.98 | 90.2 | |
5 | 10 | 9.8 | 9.8 | 100 | |
NA | 0 | 0 | NA | NA | |
total N=102 · valid N=102 · x̄=3.54 · σ=0.91 |
Da fällt eine Tabelle raus, die uns Median, die Häufigkeiten, fehlende Werte und in der Fußzeile gleich auch gesamt-N, Mittelwert und Standardabweichung anzeigt.
12.1.2 Variabilität
Hier ist die Situation ähnlich wie im vorigen Abschnitt: Varianz und Standardabweichung sind geschenkt, Variation ist seltsam.
#> [1] 0.8251796
#> [1] 0.908394
Wenig überraschend. Solange ihr euch merken könnt, dass Standardabweichung auf englisch standard deviation heißt, sollte das kein Problem sein.
Die Variation hingegen wird so selten als Maß benutzt, dass es dafür keine dedizierte Funktion gibt – analog der Modus-Sache von weiter oben. Da das allerdings wirklich kein Mensch braucht, bleibt euch nichts anderes übrig als die Varianz zu berechnen und mit \(n-1\) zu multiplizieren.
¯\_(ツ)_/¯
Okay, ihr könnt natürlich auch ganz cool sein und euch eine eigene Funktion schreiben! Wie die großen!
So etwa:
variation <- function(x) {
# NA rauskicken
x <- x[!is.na(x)]
# Variation berechnen
variation <- sum((x - mean(x))^2)
# Und zurückgeben
return(variation)
}
Das war… Definitiv etwas, das wir gerade getan haben!
#> [1] 83.34314
Funktioniert.
Braucht nur wie gesagt kein Mensch, aber jetzt habt ihr immerhin mal so grob gesehen, wie Funktionen schreiben in R funktioniert.
’Tis not witchcraft.
Aber okay, ich schweife ab.
12.1.3 Variabilität auf Ordinalskala
Auf Ordinalskala haben wir den Mittelwert nicht zur Verfügung. Das macht uns alle sehr traurig.
Aus diesem Grund müssen wir auf den Median ausweichen, den wir dann noch dezent mit Quantilen garnieren können.
Wenn wir nochmal die Variable Zufriedenheit
aus qmsurvey
nehmen:
#> [1] 4
Wenig beeindruckend.
Für Quantile:
#> 0% 25% 50% 75% 100%
#> 1 3 4 4 5
Wir können auch beliebige andere Quantile berechnen, die quantile
-Funktion gibt uns nur per default Quartile, weil die am gängigsten sind:
#> 33.3% 66.6%
#> 19 21
Über probs
geben wir der Funktion einen Vektor aus Wahrscheinlichkeiten für die entsprechenden Quantile. Hier haben wir also das 33%- und das 66%-Quantil berechnet.
Um beispielsweise Dezile zu erhalten brauchen wir also einen Vektor der Wahrscheinlichkeiten von 0 bis 10 in 10% oder 0.1er Schritten, den wir einfach via seq
erstellen können:
#> 0% 10% 20% 30% 40% 50% 60% 70% 80% 90% 100%
#> 18.0 19.0 19.0 19.0 20.0 20.0 21.0 22.0 24.8 27.0 37.0
Ein weiteres gängiges Maß für die Variabilität ist die Spannweite, die auf Ordinalskala und aufwärts sinnvoll sein kann.
Gemeint ist damit nichts anderes als der Abstand vom kleinsten zum größten Wert, also nichts besonders komplexes. In R können wir da die min
und max
-Funktionen für die beiden Werte benutzen, oder range
für beide zusammen:
#> [1] 18
#> [1] 37
#> [1] 18 37
#> [1] 19
Die Altersspanne in QM (und damit so grob im Psychologiestudiengang) ist also 19 Jahre.
Bahnbrechende Erkenntnisse die wir hier einfach so raushauen.
12.2 Bivariate Statistiken
12.2.1 Nominalskala
Auf Nominalskala ist \(\chi^2\) unser bester Freund. Leider ist es auch ein etwas übelriechender, buckliger Freund, den wir nur selten zu uns nach Hause einladen – zumindest in Deskriptivstatistikland.
Als Beispiel nehemn wir mal Geschlecht und Cannabiskonsum:
library(sjPlot)
sjt.xtab(qmsurvey$cannabis, qmsurvey$behandlung,
show.exp = TRUE, show.legend = TRUE)
Wie oft konsumierst Du Cannabis |
Warst du schon einmal in psychotherapeutischer Behandlung |
Total | |
---|---|---|---|
Ja | Nein | ||
Nie |
26 22 |
44 48 |
70 70 |
Selten |
4 8 |
22 18 |
26 26 |
Regelmäßig |
2 2 |
4 4 |
6 6 |
Total |
32 32 |
70 70 |
102 102 |
χ2=4.180 · df=2 · Cramer’s V=0.202 · Fisher’s p=0.092 |
observed values
expected values
Die Berechnung in R ist einfach:
#> Warning in chisq.test(qmsurvey$cannabis, qmsurvey$behandlung): Chi-squared
#> approximation may be incorrect
#>
#> Pearson's Chi-squared test
#>
#> data: qmsurvey$cannabis and qmsurvey$behandlung
#> X-squared = 4.1801, df = 2, p-value = 0.1237
#> Warning in chisq.test(qmsurvey$behandlung, qmsurvey$cannabis): Chi-squared
#> approximation may be incorrect
#>
#> Pearson's Chi-squared test
#>
#> data: qmsurvey$behandlung and qmsurvey$cannabis
#> X-squared = 4.1801, df = 2, p-value = 0.1237
Das Output ist nicht so hübsch und etwas clunky und beinhaltet auch p-Wert und Freiheitsgrade, die euch eigentlich nur in einem inferenzstatistischen Kontext interessieren.
Darüber hinaus kann es sein, dass der berechnete Wert nicht mit händischer Berechnung übereinstimmt.
Der Grund dafür ist Yates’ Kontinuitätskorrektor, die nur bei 2x2-Tabellen greift. Solltet ihr die Korrektur explizit deaktivieren wollen, könnt ihr das Argument correct = FALSE
setzen.
Für kompakteres Output haben wir eine entsprechende Funktion in der tadaatoolbox:
#> [1] 4.180078
Als nächstes haben wir die auf \(\chi^2\)-basierenden Statistiken Phi (\(\phi\)) und Cramer’s V, die entsprechenden Funktionen kommen auch aus der tadaatoolbox:
#> [1] 0.202438
#> [1] 0.202438
Oh nein, Phi funktioniert nicht. Wie unerwartet.
Ihr mögt euch erinnern, dass Phi nur für 2x2-Tabellen definiert ist, wobei Cramer’s V bei beliebigen Tabellendimensionen funktioniert – darüberhinaus sind die beiden Statistiken im Falle von 2x2-Tabellen auch noch identisch, weshalb sich die Frage aufwirft, wieso überhaupt jemand Phi benutzen sollte – aber nun gut, die Funktionen sind da.
Als letztes gibt es noch den Kontingenzkoeffizient C – keine Ahnung was der für eine Daseinsberechtigung hat, aber wenn ihr ihn mal braucht:
#> [1] 0.1984133
Sagte ich als letztes? Das war gelogen. Es gibt noch Lambda (\(\lambda\)) in drei Geschmacksrichtungen:
- “Sorum” (Variablen wie definiert)
- “Andersrum” (x- und y-Variable vertauschen)
- Symmetrisch
#> [1] 0
#> [1] 0
#> [1] 0
Es gibt zusätzlich die Funktion assocstats
aus dem vcd
-Package. Da gebt ihr allerdings nicht beide Variablen einzeln an, sondern einen table
der beiden Variablen:
#>
#> Nie Selten Regelmäßig
#> Ja 26 4 2
#> Nein 44 22 4
#> X^2 df P(> X^2)
#> Likelihood Ratio 4.5754 2 0.10150
#> Pearson 4.1801 2 0.12368
#>
#> Phi-Coefficient : NA
#> Contingency Coeff.: 0.198
#> Cramer's V : 0.202
So könnt ihr die gängigen Statistiken auf einen Blick sehen.
Alternativ haben wir dafür natürlich auch was in der tadaatoolbox:
\(\chi^2\) | Cramer’s V | \(\lambda_x\) | \(\lambda_y\) | \(\lambda_{xy}\) | c |
---|---|---|---|---|---|
4.18 | 0.2 | 0 | 0 | 0 | 0.2 |
Das Argument print = "markdown"
ist nur für das Output in diesem Buch bzw. in RMarkdown-Dokumenten relevant, für die Nutzung in der Konsole könnt ihr das weglassen.
Aber ihr benutzt ja jetzt eh alle brav RMarkdown, seit ihr das entsprechende Kapitel dazu hier überflogen habt, ja?
Ja.
Brav.
12.2.2 Ordinalskala
Auf Ordinalskala haben wir primär Gamma, Somers’ D und eine Handvoll Tau abzufrühstücken.
Um’s einach zu halten bedienen wir uns auch hier wieder der tadaatoolbox:
#> [1] 0.1170716
#> [1] 0.1090856
#> [1] 0.07730777
#> [1] 0.09319667
Wobei Somers’ D analog Lambda drei Geschmacksrichtungen hat.
Oder auch hier, alle zusammen:
\(\gamma\) | \(D_x\) | \(D_y\) | \(D_{xy}\) | \(\tau_A\) | \(\tau_B\) | \(\tau_C\) |
---|---|---|---|---|---|---|
0.12 | 0.11 | 0.08 | 0.09 | 0.07 | 0.09 | 0.09 |
Da sind dann auch gleich die Taus mit drin. Die könnt ihr auch einzeln bekommen:
#> [1] 0.07202485
#> [1] 0.09183225
#> [1] 0.0891484
Da \(\tau_b\) auch als Korrelationskoeffizient durchgeht, gibts’s den auch via cor
:
#> [1] 0.09183225
12.2.3 Intervallskala
Hier fangen die spaßigen Dinge an. Korrelationen und so.
Und was ist Korrelation?
Die standardisierte Kovarianz, nämlich.
Und was ist Kovarianz?
#> [1] 2.810619
Nämlich.
cov
gibt uns die Kovarianz, wenig überraschend, und über use = "complete.obs"
spezifizieren wir, dass Beobachtungspaare mit fehlenden Werten (NA
) ausgelassen werden sollen. Eine Alternative wäre "pairwise.complete.obs"
, was in diesem Fall keinen Unterschied macht, aber prinzipiell wollt ihr eins von beidem benutzen um zu verhindern, dass ihr als Ergebnis nur ein trauriges NA
bekommt.
Als nächstes also die Korrelationskoeffizienten:
- Pearson’s \(r\)
- Spearman’s \(\rho\)
Beide bekommen wir mit cor
, abhängig vom method
-Argument. Der Default ist "pearson"
#> [1] 0.5225806
#> [1] 0.5303142
Über method = "kendall"
bekommen wir auch Kendall’s \(\tau_b\). Just in case.
12.2.4 Weitere Statistiken
Es gibt noch viele andere Statistiken, die für euch mehr oder weniger spannend sind. Viele andere sind für euch sogar vollkommen irrelevant.
Point being: Das hier
12.2.4.1 \(\eta^2\)
Ihr lernt \(\eta^2\) in QM1 als deskriptive Statistik kennen, wo es als nicht-lineare Alternative zu \(r\) oder \(\rho\) verwendet wird. Das ist auch okay so, aber generell findet \(\eta^2\) eher als Effektgröße bei der ANOVA Anwendung, weshalb auch R eher diesen Kontext erwartet.
Wenn ihr \(\eta^2\) aber “einfach so” haben wollt, braucht ihr etwas Umweg:
#> term etasq
#> 1 alkohol 0.09
Wir haben hier zuerst ein ANOVA-Modell mittels aov
erstellt, wo wir die Variable alter
aufgeteilt haben in Gruppen abhängig von der Variable alkohol
. Dann haben wir eta_sq
aus sjstats
für die Berechnung von \(\eta^2\) benutzt, und… sind damit fertig.
Als Bonus-Feature haben wir da seit der Version 0.16.0
auch was passendes in der tadaatoolbox
:
#> [1] 0.08971183
Unter der Haube benutzen wir da übrigens DescTools::EtaSq
. Ihr seht, es gibt da etliche Möglichkeiten, wie so oft in R.
12.2.4.2 Cohen’s Kappa (\(\kappa\))
#> rater1 rater2 rater3
#> 1 3 3 2
#> 2 3 6 1
#> 3 3 4 4
#> 4 4 6 4
#> 5 5 2 3
#> 6 5 4 2
#> Cohen's Kappa for 2 Raters (Weights: unweighted)
#>
#> Subjects = 20
#> Raters = 2
#> Kappa = 0.119
#>
#> z = 1.16
#> p-value = 0.245
12.2.4.3 Kendall’s \(W\)
#> Kendall's coefficient of concordance W
#>
#> Subjects = 20
#> Raters = 3
#> W = 0.502
#>
#> Chisq(19) = 28.6
#> p-value = 0.0724
#>
#> Coefficient may be incorrect due to ties