Nachdem man eine Datei korrekt importiert hat (d.h. wenn die Daten in einem Data Frame bzw. Tibble vorhanden sind und alle Spalten den passenden Datentyp haben), kann man mit der statistischen Analyse beginnen. Der erste Schritt dabei ist meist, sich mit deskriptiven Statistiken eine Beschreibung der vorhandenen Daten zu erzeugen.
Als Beispiel dazu importieren wir die Datei lecturer.dat aus der letzten Übung:
Die name-Spalte brauchen wir für unsere nachfolgenden Betrachtungen nicht mehr, deswegen entfernen wir sie:
df$name =NULL
Es gibt in R eine Reihe an Funktionen, welche zusammenfassende Statistiken einer Variablen (bzw. eines Vektors) berechnen. Nützliche Funktionen sind z.B. mean(), sd(), var(), min(), max(), median(), range() und quantile(). Den Mittelwert eines Vektors (und damit auch einer Spalte von df) kann man also wie folgt berechnen:
mean(df$friends)
[1] 7.8
Vor allem bei Faktoren ist es interessant zu wissen, wie viele Stufen (unterschiedliche Werte) diese beinhalten. Auch bei numerischen oder Text-Vektoren mit mehrfach vorkommenden Werten ist es manchmal praktisch, die unterschiedlichen Werte herauszufinden. Dafür gibt es die Funktion unique(), welche die einzigartigen Elemente eines Vektors zurückgibt:
unique(c(1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2))
[1] 1 2
unique(df$friends)
[1] 5 2 0 4 1 10 12 15 17
unique(df$job)
[1] Lecturer Student
Levels: Lecturer Student
All diese Berechnungen müsste man nun für jede interessierende Spalte wiederholen (weil diese nur mit Vektoren funktionieren), was relativ mühsam wäre. Deswegen gibt es die Funktion sapply(), welche eine Funktion auf jede Spalte eines Data Frames einzeln anwendet. Möchte man also den Mittelwert für jede numerische Spalte von df berechnen, kann man dies so tun:
sapply(df[, -c(1, 2)], mean)
friends alcohol income neurotic
7.8 17.6 18511.0 13.1
Tipp
Für die Berechnung der Spaltenmittelwerte kann man auch die uns bereits bekannte Funktion colMeans() verwenden.
Beachten Sie, dass df[, -c(1, 2)] alle Spalten aus df ausgenommen der ersten und zweiten bedeutet. So kann man jede beliebige Funktion auf einzelne Spalten anwenden, z.B. die Standardabweichung:
sapply(df[, -c(1, 2)], sd)
friends alcohol income neurotic
6.142746 7.042727 18001.354548 3.984693
Es gibt aber auch spezielle Funktionen, welche direkt mehrere statistische Kenngrößen für alle Spalten eines Data Frames berechnen. Im Folgenden gehen wir näher auf drei dieser Funktionen ein, nämlich summary(), describe() und stat.desc(). Nur summary() ist standardmäßig bei R dabei, die anderen zwei Funktionen erfordern die Installation von zusätzlichen Paketen.
Die Funktion summary()
Die Funktion summary() liefert eine geeignete Zusammenfassung für jede Spalte eines Data Frames (Tibbles). Numerische Spalten sowie Datumsspalten werden mit sechs Werten beschrieben:
Minimum
Erstes Quartil (25%-Perzentil)
Median (50%-Perzentil)
Drittes Quartil (75%-Perzentil)
Maximum
Arithmetischer Mittelwert
Für Faktoren werden die Stufen sowie die Anzahl an Fällen pro Stufe aufgelistet.
summary(df)
birth_date job friends alcohol income neurotic
Min. :1949-10-10 Lecturer:5 Min. : 0.0 Min. : 5.00 Min. : 10 Min. : 7.00
1st Qu.:1971-04-01 Student :5 1st Qu.: 2.5 1st Qu.:15.25 1st Qu.: 3500 1st Qu.:10.75
Median :1975-06-27 Median : 7.5 Median :17.50 Median :15000 Median :13.00
Mean :1975-12-17 Mean : 7.8 Mean :17.60 Mean :18511 Mean :13.10
3rd Qu.:1984-08-10 3rd Qu.:12.0 3rd Qu.:20.00 3rd Qu.:31750 3rd Qu.:14.00
Max. :1989-01-23 Max. :17.0 Max. :30.00 Max. :50000 Max. :21.00
Auch bei anderen Objekten als Data Frames liefert summary() oft eine sinnvolle kurze Beschreibung.
Die Funktion describe()
Eine weitere Möglichkeit, noch mehr statistische Kenngrößen für numerische Spalten auszugeben, bietet die Funktion describe() aus dem psych-Paket. Nicht-numerische Spalten werden hier nicht vernünftig zusammengefasst, deshalb sollte man der Funktion nur numerische Spalten übergeben.
library(psych)describe(df[, 3:6])
vars n mean sd median trimmed mad min max range skew kurtosis se
friends 1 10 7.8 6.14 7.5 7.62 7.41 0 17 17 0.11 -1.75 1.94
alcohol 2 10 17.6 7.04 17.5 17.62 3.71 5 30 25 -0.03 -0.75 2.23
income 3 10 18511.0 18001.35 15000.0 16887.50 19940.97 10 50000 49990 0.45 -1.48 5692.53
neurotic 4 10 13.1 3.98 13.0 12.88 2.97 7 21 14 0.36 -0.66 1.26
Hinweis
Eigentlich kann man Funktionen aus Paketen auch ohne Aktivierung verwenden. Dazu muss man den Paketnamen gefolgt von :: voranstellen. Im vorigen Beispiel könnte man also auf library(psych) verzichten und die Funktion mit psych::describe() trotzdem verwenden.
Man kann diese Funktion auch auf einzelne Gruppen separat anwenden. Im Beispiel könnte man dies getrennt für jede Stufe von df$job tun. Dazu verwendet man die verwandte Funktion describeBy():
Descriptive statistics by group
group: Lecturer
vars n mean sd median trimmed mad min max range skew kurtosis se
friends 1 5 2.4 2.07 2 2.4 2.97 0 5 5 0.11 -2.03 0.93
alcohol 2 5 16.0 9.62 15 16.0 7.41 5 30 25 0.28 -1.72 4.30
income 3 5 33400.0 12561.85 35000 33400.0 19273.80 20000 50000 30000 0.10 -1.98 5617.83
neurotic 4 5 15.0 4.18 14 15.0 4.45 10 21 11 0.25 -1.72 1.87
------------------------------------------------------------------------------------------
group: Student
vars n mean sd median trimmed mad min max range skew kurtosis se
friends 1 5 13.2 2.77 12 13.2 2.97 10 17 7 0.23 -1.89 1.24
alcohol 2 5 19.2 3.56 18 19.2 2.97 16 25 9 0.66 -1.43 1.59
income 3 5 3622.0 4135.69 3000 3622.0 4299.54 10 10000 9990 0.48 -1.64 1849.54
neurotic 4 5 11.2 3.03 13 11.2 1.48 7 14 7 -0.37 -2.01 1.36
Das erste Argument ist das zu beschreibende Data Frame, und das zweite Argument ist die Spalte, nach der gruppiert werden soll (üblicherweise ein Faktor).
Die Funktion stat.desc()
Das Paket pastecs beinhaltet die Funktion stat.desc() zur Beschreibung von Daten. Mit der Funktion round() sollte man hier außerdem einstellen, wie viele Kommastellen ausgegeben werden sollen, da die Ausgabe der Funktion sonst relativ unübersichtlich ist. Wenn das Argument norm=TRUE gesetzt wird, werden für alle Spalten Tests auf Normalverteilung durchgeführt.
In R kann man Dezimalzahlen auch in der sogenannten wissenschaftlichen Notation anschreiben. Hier verwendet man eine Darstellung mit Zehnerpotenzen, die man mit e eingeben kann – e kann man als “mal zehn hoch” lesen.
1e0# 1 mal 10 hoch 0
[1] 1
-4e0# -4 mal 10 hoch 0
[1] -4
1e1# 1 mal 10 hoch 1
[1] 10
3.5e2# 3.5 mal 10 hoch 2
[1] 350
1e-2# 1 mal 10 hoch -2
[1] 0.01
1e-15# 1 mal 10 hoch -15 = 0.000000000000001
[1] 1e-15
Tipp
Möchte man die Ausgabe von Zahlen in der wissenschaftlichen Notation weitgehend vermeiden, kann man zu Beginn einer Sitzung folgende Option setzen:
options(scipen=100)
Eine andere Möglichkeit ist es, Zahlen mit format() als normale Dezimalzahlen darzustellen. Dies ergibt aber immer einen Character-Vektor, erfordert dafür aber nicht das Ändern einer Option:
format(1e-15, scientific=FALSE)
[1] "0.000000000000001"
Gruppieren mit by()
Für die Funktion stat.desc() gibt es keine direkte Variante für gruppierte Daten (so wie describeBy() für describe()). Es gibt aber in R die generische Funktion by(), welche beliebige Funktionen auf gruppierte Daten anwendet. Das erste Argument ist hier wie üblich der Datensatz, das zweite Argument ist die Gruppierungsspalte, und das dritte Argument ist die Funktion, die auf die gruppierten Daten angewendet werden soll.
df$job: Lecturer
vars n mean sd median trimmed mad min max range skew kurtosis se
friends 1 5 2.4 2.07 2 2.4 2.97 0 5 5 0.11 -2.03 0.93
alcohol 2 5 16.0 9.62 15 16.0 7.41 5 30 25 0.28 -1.72 4.30
income 3 5 33400.0 12561.85 35000 33400.0 19273.80 20000 50000 30000 0.10 -1.98 5617.83
neurotic 4 5 15.0 4.18 14 15.0 4.45 10 21 11 0.25 -1.72 1.87
------------------------------------------------------------------------------------------
df$job: Student
vars n mean sd median trimmed mad min max range skew kurtosis se
friends 1 5 13.2 2.77 12 13.2 2.97 10 17 7 0.23 -1.89 1.24
alcohol 2 5 19.2 3.56 18 19.2 2.97 16 25 9 0.66 -1.43 1.59
income 3 5 3622.0 4135.69 3000 3622.0 4299.54 10 10000 9990 0.48 -1.64 1849.54
neurotic 4 5 11.2 3.03 13 11.2 1.48 7 14 7 -0.37 -2.01 1.36
Möchte man der Funktion im dritten Argument selbst Argumente übergeben (z.B. norm=TRUE für die Funktion stat.desc()), kann man dies mit weiteren Argumenten ganz am Ende tun:
df$job: Lecturer
friends alcohol income neurotic
Min. :0.0 Min. : 5 Min. :20000 Min. :10
1st Qu.:1.0 1st Qu.:10 1st Qu.:22000 1st Qu.:13
Median :2.0 Median :15 Median :35000 Median :14
Mean :2.4 Mean :16 Mean :33400 Mean :15
3rd Qu.:4.0 3rd Qu.:20 3rd Qu.:40000 3rd Qu.:17
Max. :5.0 Max. :30 Max. :50000 Max. :21
------------------------------------------------------------------------------------------
df$job: Student
friends alcohol income neurotic
Min. :10.0 Min. :16.0 Min. : 10 Min. : 7.0
1st Qu.:12.0 1st Qu.:17.0 1st Qu.: 100 1st Qu.: 9.0
Median :12.0 Median :18.0 Median : 3000 Median :13.0
Mean :13.2 Mean :19.2 Mean : 3622 Mean :11.2
3rd Qu.:15.0 3rd Qu.:20.0 3rd Qu.: 5000 3rd Qu.:13.0
Max. :17.0 Max. :25.0 Max. :10000 Max. :14.0
Die Funktion stat.desc() liefert bereits das Ergebnis des Shapiro-Wilk-Tests auf Normalverteilung (die Einträge normtest.W und normtest.p enthalten den Wert der Teststatistik bzw. die Signifikanz). Wenn normtest.p signifikant ist (z.B. kleiner als 0.05), dann kann man die Nullhypothese der Normalverteilung verwerfen. Man kann den Shapiro-Wilk-Test auch direkt mit der Funktion shapiro.test() aufrufen:
shapiro.test(df$income)
Shapiro-Wilk normality test
data: df$income
W = 0.89721, p-value = 0.2041
Mit der by()-Funktion kann man den Test auch getrennt auf verschiedene Gruppen anwenden.
by(df$income, df$job, shapiro.test)
df$job: Lecturer
Shapiro-Wilk normality test
data: dd[x, ]
W = 0.93425, p-value = 0.6256
------------------------------------------------------------------------------------------
df$job: Student
Shapiro-Wilk normality test
data: dd[x, ]
W = 0.89439, p-value = 0.3796
Der Kolmogorov-Smirnov-Test kann gegebene Daten auf beliebige Verteilungen testen, d.h. natürlich auch auf Normalverteilung. Im Falle der Normalverteilung ist aber der Shapiro-Wilk-Test vorzuziehen, da dieser speziell auf die Normalverteilung zugeschnitten ist und daher mehr statistische Power besitzt.
Exact one-sample Kolmogorov-Smirnov test
data: df$income
D = 0.18182, p-value = 0.8389
alternative hypothesis: two-sided
Da die Stichprobengröße in unserem Beispiel nur sehr klein ist, lassen sich aber ohnehin keine vernünftigen Aussagen über die Verteilung der Daten treffen.
Test auf Varianzhomogenität
Der Levene-Test prüft auf Gleichheit der Varianzen (Homoskedastizität) von zwei oder mehr Gruppen. Die Nullhypothese ist, dass die Varianzen in allen Gruppen gleich sind. In R führt man den Test mit der Funktion leveneTest() aus dem Paket car durch. Dazu sehen wir uns die Beispieldaten Moore an, welche mit dem Paket car automatisch geladen werden.
partner.status conformity fcategory fscore
42 high 13 high 57
43 high 16 low 35
44 high 10 high 52
45 high 15 medium 44
summary(Moore)
partner.status conformity fcategory fscore
high:23 Min. : 4.00 high :15 Min. :15.00
low :22 1st Qu.: 8.00 low :15 1st Qu.:35.00
Median :12.00 medium:15 Median :43.00
Mean :12.13 Mean :43.11
3rd Qu.:15.00 3rd Qu.:55.00
Max. :24.00 Max. :68.00
Der Levene-Test für die Spalte conformity gruppiert nach der Spalte fcategory wird wie folgt aufgerufen:
leveneTest(Moore$conformity, Moore$fcategory)
Levene's Test for Homogeneity of Variance (center = median)
Df F value Pr(>F)
group 2 0.046 0.9551
42
In diesem Beispiel kann die Nullhypothese der Varianzgleichheit der Variable conformity über die Gruppen in fcategory also nicht verworfen werden, da der p-Wert (hier als Pr(>F) bezeichnet) mit 0.9551 extrem groß ist.
Übungen
Übung 1
Berechnen Sie statistische Kenngrößen wie Mittelwert, Median, Minimum und Maximum für die vier numerischen Variablen Global_active_power, Global_reactive_power, Voltage und Global_intensity aus den Daten die Sie in der vorigen Einheit bereits importiert haben (Individual Household Electric Power Consumption).
Berechnen Sie die Kenngrößen mit der Funktion sapply().
Berechnen Sie die Kenngrößen mit der Funktion summary().
Berechnen Sie die Kenngrößen mit der Funktion describe() aus dem Paket psych.
Wenden Sie die Funktion stat.desc() aus dem Paket pastecs auf die Daten an (runden Sie die Ergebnisse auf eine Nachkommastelle).
Wie groß ist die gemessene mittlere Spannung Voltage?
Wie groß ist der Median der globalen Wirkleistung Global_active_power?
Wie viele fehlende Werte gibt es insgesamt bzw. pro Spalte?
Übung 2
Das Paket palmerpenguins beinhaltet im Tibble penguins Messdaten über physische Eigenschaften von drei Pinguinarten (Adélie, Chinstrap und Gentoo). Aktivieren Sie das Paket und beantworten Sie mit geeigneten R-Befehlen folgende Fragen:
Aus wie vielen Zeilen und Spalten besteht der Datensatz penguins?
Die Spalten species, island und sex sind Faktoren – wie viele Stufen gibt es jeweils?
Berechnen Sie deskriptive Statistiken getrennt für die drei Pinguinspezies! Wie lauten die Mittelwerte der Spalten bill_length_mm, bill_depth_mm und flipper_length_mm für die drei Spezies?
Gibt es fehlende Werte in den Daten?
Übung 3
Im Datensatz mtcars befinden sich diverse Kenngrößen von 32 Autos.
Berechnen Sie das Minimum und Maximum sowie den Mittelwert und den Median von allen Spalten.
Überprüfen Sie, ob die Daten in der Spalte mpg normalverteilt sind.
Überprüfen Sie, ob Varianzhomogenität für die Daten in der Spalte mpg gegeben ist für die drei Gruppen, die durch die Spalte cyl definiert sind.
Übung 4
Laden Sie die Beispieldaten lecturer.dat und berechnen Sie die Mittelwerte aller numerischen Spalten gruppiert nach der Spalte job. Verwenden Sie dazu die Funktion by(). Warum funktioniert diese Berechnung nicht mit der Funktion mean?
Übung 5
Erzeugen Sie mit der Funktion rnorm() einen Vektor mit 5001 standardnormalverteilten Zufallszahlen und überprüfen Sie, ob dieser Vektor normalverteilt ist. Verwenden Sie dazu den Shapiro-Wilk-Test. Wie interpretieren Sie das Ergebnis?