---
title: "L2: Abbildungen in R"
author: "Jonathan Harrington / Ulrich Reubold"
date: "30. April / 07. Mai 2019"
output: html_document
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```

Bibliotheken und Data-Frames laden:
```{r,eval=FALSE,message=FALSE}
# ggf.: 
install.packages('ggplot2')
install.packages('gridExtra')
```
```{r,echo=FALSE}
pfad = "/Users/reubold/Documents/R"
```

```{r,message=FALSE}
pfadu = "http://www.phonetik.uni-muenchen.de/~jmh/lehre/Rdf"
library(ggplot2)
#VOR der Vorlesung:
#Alles zwischen ### ### an einem Ort mit stabiler Internetverbindung durchführen

###
asp = read.table(file.path(pfadu, "asp.txt"))
coronal = read.table(file.path(pfadu, "coronal.txt"))
int.df = read.table(file.path(pfadu, "intdauer.txt"))
v.df = read.table(file.path(pfadu, "vdata.txt"))
save(list=c("asp","coronal","int.df","v.df"),file=file.path(pfad, "L2.RData"))
###

#dann
load(file.path(pfad, "L2.RData"))
```

##0. Einleitung
Zur Erinnerung, `df` ist ein Data-Frame.

Ein Data-Frames besteht aus 

* Reihen (Beobachtungen) 

und

* Spalten (Variablen)

```{r}
# Prüfen, ob es ein Data-Frame ist:
class(asp)
# die ersten paar Beobachtungen:
head(coronal)
```

#### Zugriff auf Elemente
```{r,eval=FALSE}
coronal[r,s]  #Zugriff auf Reihe(n) r und Spalte(n) s
```
also z.B.:
```{r}
coronal[2,2:4]
```

## 1. Numerische und kategoriale Variablen
Die Spalten können **numerische** und **kategoriale** Variablen enthalten. 

#### Numerische Variablen: Ganzzahlen (`integer`) oder kontinuierlich (`numeric`: Gleitkommazahlen):

```{r}
asp$d[1:5]
class(asp$d)
# oder
with(asp, class(d))

int.df$Dauer[1:5]
class(int.df$Dauer)
```

#### Kategoriale Variablen sind Faktoren mit Stufen:
```{r}
class(coronal$Socialclass)
# [1] "factor"

# Die Stufen davon. z.B. die ersten paar anschauen...
coronal$Socialclass[1:10]
# oder
levels(coronal$Socialclass)
```

Die Inhalte von Faktoren können Buchstaben oder Zahlen sein; diese werden aber grundsätzlich anders behandelt als numerische Daten oder Zeichenketten, da bei der Klasse ` factor ` wiederholt auftauchende Symbole als einer Kategorie angehörig erkannt werden.

## 2. Typische Untersuchungen in der Phonetik

### Wird y von x beeinflusst?

#### 1. y = numerisch, x = kategorial

* 1.1 z.B. *unterscheiden sich /i, e, a/ in der Dauer?*
* 1.2 = **Inwiefern wird die Dauer (y) vom Vokal (x) beeinflusst?**
* 1.3 Abbildungen dafür: `geom_boxplot()` 
    + auch: `geom_histogram()`  und `stat_density()`

#### 2. y = kategorial, x = kategorial

* 2.1 z.B. *Wörter wie Sohn, Sonne... können mit /s/ oder /z/ produziert werden. Kommt /s/ eher in Bayern und /z/ eher in Hamburg vor?*
* 2.2 = **Inwiefern wird die Stimmhaftigkeit (y: ob /s/ oder /z/) vom Dialekt (x: oben Bayern oder Hamburg) V) beeinflusst?**
* 2.3 Abbildungen dafür: `geom_barchart()`

#### 3. y = numerisch, x = numerisch 

* 3.1 z.B. *gibt es eine größere Dauer bei einer größeren Mundöffnung?*
* 3.2 **Inwiefern wird die Dauer (y) von der Mundöffnung (x) beeinflusst?**
* 3.3 Abbildung dafür: xy-Plot über `geom_point()` und/oder `geom_line()`

### Zusammenfassung: 

* 1. y = numerisch, x = kategorial: `geom_boxplot()`.  Auch `geom_histogram()`, `stat_density()`

* 2. y = kategorial, x = kategorial: `geom_bar()`

* 3. y = numerisch, x = numerisch: `geom_point()`, `geom_line()`

#### Siehe auch die Handbücher zu 

* R und Abbildungen im Allgemeinen http://www.cookbook-r.com/Graphs/
* zu ggplot2 im Speziellen (gg = grammar of graphics) http://docs.ggplot2.org/current/

Syntax in ggplot2:
```{r,eval=FALSE}
A + 
  B + 
  C + 
  D + 
  ...
```

`A`, `B`, `C`... sind modulare Befehle, die durch `+` aneinandergehängt werden.

Hier:
`A. (data-frame) + B. (Variablen) + C. (Art der Abbildung)`

Die oben genannten Methoden `geom_boxplot()`, `geom_bar()`, `geom_point()` usw. sind Arten der Abbildung.

Was in dieser Form noch abschreckend wirken mag, wird durch konkrete Anwendungen in den folgenden Beispielen schnell deutlicher.

## 3. `geom_boxplot()`: y = numerisch, x = kategorial
```{r}
head(asp)
```
Inwiefern wird die Aspirationsdauer (d) von der Artikualtionsstelle (Kons) beeinflusst?

* y: d    (numerisch)
* x: Kons (kategorial)

```{r}
ggplot(asp) + 
  aes(y = d, x = Kons) + 
  geom_boxplot()
```
```{r,eval=FALSE}
# oder in getrennten Komponenten ablegen.
# A
p1 = ggplot(asp)
# B
p2 = aes(y = d, x = Kons)
# C
p3 = geom_boxplot()
# A + B + C
p1 + p2 + p3
# oder A + B + C ablegen
erg = p1 + p2 + p3
# Bild
erg
```

Weitere Informationen zum Boxplot:
Die Linie ist der **Median**;  der 'Box' ist der **interquartale Bereich**.

### Exkurs: Median und Interquartaler Bereich

* *Median* (oder *50% Quantil*): **Der mittlere Wert in der sortierten Reihenfolge**:

```{r}
x = c(15, 19, -1, 10, 11, 18, 90000)
median(x)

# das gleiche
quantile(x, .5)

# Nach Sortierung sieht man, dass 15 an der mittleren Stelle vorkommt.
sort(x)
```

N.B.: Wenn es eine gerade Anzahl an Beobachtungen gibt, ist das arithmetische Mittel der beiden mittleren Zahlen aus einem sortierten Zahlenvektor der Median-Wert!

```{r}
x2 = c(15, 19, 13,-1, 10, 11, 18, 90000)
median(x2)
sort(x2)
mean(c(13,15))
```


Es gibt auch andere Bruchteile von Quantilen. Insbesondere wird von dem interquartilen Bereich Gebrauch gemacht, um die Streuung einzuschätzen:

```{r}
IQR(x)

# das Gleiche
quantile(x, .75) - quantile(x, .25)
#da
quantile(x, .75)
quantile(x, .25)
```


## 4. `geom_bar()`: y ist kategorial, x ist kategorial
```{r}
head(coronal)
```
Inwiefern wird die Artikulationsstelle (Fr) von der Region (Region) beeinflusst?

* y: Fr (kategorial)
* x: Region (kategorial)

```{r}
p1 = ggplot(coronal)
p2 = aes(fill = Fr, x = Region)
# Um Häufigkeiten darzustellen
p3 = geom_bar()
p1 + p2 + p3
# Um die Bars der Stufen von x nebeneinander zu stellen
p4 = geom_bar(position="dodge")
p1 + p2 + p4

# Um Proportionen darzustellen
p5 = geom_bar(position="fill")
p1 + p2 + p5
```
```{r,eval=FALSE}
#oder
ggplot(coronal) +
  aes(fill = Fr, x = Region) +
  geom_bar(position="fill")
```





## 5. `geom_point()`, `geom_line()`:  y ist numerisch, x ist numerisch
```{r}
head(int.df)
```
Inwiefern wird die Dauer (Dauer) von der Intensität (dB) beeinflusst?

* y: Dauer (numerisch) 
* x: dB (numerisch)

```{r,eval=FALSE}
# Nur Linie
ggplot(int.df) +  
  aes(x = dB, y = Dauer) + 
  geom_line() 
```
```{r}
#oder: der aes()-Teil kann auch direkt bei ggplot() eingegeben werden:
ggplot(int.df,aes(x = dB, y = Dauer)) +  
  geom_line() 
# Nur Punkte
ggplot(int.df, aes(x = dB, y = Dauer)) + 
  geom_point() 
# Beide
ggplot(int.df, aes(x = dB, y = Dauer)) + 
  geom_line() + 
  geom_point()
```

## 6. Achsenbeschriftung: `+ xlab()` `+ ylab()` `+ ggtitle()`
```{r}
# Früherer Boxplot-Befehl
p1 = ggplot(asp) + 
  aes(y = d, x = Kons) + 
  geom_boxplot()
# Etikettierung für die x-Achse
p2 = xlab("Artikulationsstelle")
# Etikettierung für die y-Achse
p3 = ylab("Dauer (ms)")
# Titel
p4 = ggtitle("Boxplot-Daten")
p1 + p2 + p3 + p4

# Früherer Barchart-Befehl
bar.p = ggplot(coronal) + 
  aes(x = Region, fill = Fr) + 
  geom_bar(position = "fill")
x.p = xlab("Region")
y.p = ylab("Proportion")
t.p = ggtitle("Proportionale Aufteilung von Frikativen")
bar.p + x.p + y.p + t.p
```

## 7. Bereiche: `+coord_cartesian(xlim=c(,),ylim=c(,))`

```{r}
# Früherer geom_bar() Befehl
p1 = ggplot(int.df, aes(dB, Dauer)) + geom_point() 
# xlim und ylim
p2 = coord_cartesian(xlim=c(10,60),ylim=c(30,280))
p1 + p2
```

Vorsicht: man könnte hier auch einfach `+xlim(10,60) + ylim(30,280)` schreiben, aber auf diese Weise werden die nicht-geplotteten Daten einfach ignoriert, was problematisch ist, wenn man Funktionen auf die Daten anwenden will, z.B. Regressionslinien überlagern (siehe Vorlesung Ende Mai);
siehe z.B.: https://i.stack.imgur.com/ciyvf.png

## 8. Farben
Farbenauswahl: http://www.stat.columbia.edu/~tzheng/files/Rcolor.pdf, auch:
```{r,eval=FALSE}
colors()
```

#### `geom_boxplot()`
```{r}
########## Der vorige Befehl
ggplot(asp) + 
  aes(y = d, x = Kons) + 
  geom_boxplot()
# Default-Farben
# Jeder Box mit unterschiedlichen Farben gefüllt
ggplot(asp) + 
  aes(y = d, x = Kons, fill = Kons) + 
  geom_boxplot()
# Jeder Box mit verschiedenen Umrissfarben
ggplot(asp) + 
  aes(y = d, x = Kons, col = Kons) + 
  geom_boxplot()

# Eigene Farben wählen
farben = c("green", "red")
# gefüllt
ggplot(asp) + aes(y = d, x = Kons) + geom_boxplot(fill = farben)
# Umriss
ggplot(asp) + aes(y = d, x = Kons) + geom_boxplot(col = farben)
```

Merke: 

* `col`=Umrissfarbe
* `fill`=Füllungsfarbe

#### `geom_bar()`
```{r}
########## Der vorige Befehl
p1 = ggplot(coronal) + 
  aes(x = Region, fill = Fr) + 
  geom_bar()
p1
# Eigene Farben wählen
farben = c("yellow", "green")
p2 = scale_fill_manual(values = farben) 
p1 + p2
```

N.B.: Bei `geom_boxplot()` und `geom_point()` bzw. `geom_line()` verwendet man anstelle von `scale_fill_manual(values = farben)` `scale_color_manual(values = farben)` (wenn denn in `aes(...,...,col = X2)` vorkam).

#### `geom_point()`, `geom_line()`

* `col`: Farbe. 
* `pch`: plotting character. 
Siehe http://www.endmemo.com/program/R/pchsymbols.php
* `cex`: character expansion: `cex = 2` bedeutet zwei Mal die übliche Größe
* `lwd`: Linienstärke

```{r}
########## Der vorige Befehl
ggplot(int.df, aes(x = dB, y = Dauer)) +  
  geom_point() + 
  geom_line()

#eigene Farben wählen:
ggplot(int.df, aes(x = dB, y = Dauer)) +  
  geom_point(col="purple", pch=0, cex=2) + 
  geom_line(col = "pink")
# lwd: Linienstärke
ggplot(int.df, aes(x = dB, y = Dauer)) +  
  geom_point(col="purple", pch=0, cex=2) + 
  geom_line(col = "pink", lwd=2)
```

## 9. Schriftzeichengröße
Der Default-Schriftzeichen-Größe der Achsenbeschriftung und Titel ist 11. Für die Legende etwas weniger (eventuell 10).

```{r}
p1 = ggplot(asp) + 
  aes(y = d, x = Kons) + 
  geom_boxplot() + 
  xlab("Artikulationsstelle") + 
  ylab("Dauer (ms)") + 
  ggtitle("Boxplot-Daten")
p1

# Keine Änderung
p11 = theme(text = element_text(size=11))
p1 + p11

# Schriftzeichen-Größe 16
p16 = theme(text = element_text(size=16))
p1 + p16

# Nur Achsenbeschriftung vergrößern
q24 = theme(axis.text = element_text(size=24))
p1 + q24

# Achsenbeschriftung und Bildetikettierungen unterschiedlich vergrößern
p30 = theme(text = element_text(size=30))
p1 + q24 + p30
```

## 10. Mehrere Bilder nebeneinander

Hierfür muss ein anderes Package geladen werden; es gibt allerdings [mehrere Packages für diesen Zweck](http://www.sthda.com/english/articles/24-ggpubr-publication-ready-plots/81-ggplot2-easy-way-to-mix-multiple-graphs-on-the-same-page/), z.B. `gridExtra`:

```{r}
library(gridExtra)

pic1 = ggplot(asp, aes(y = d, x = Kons))  + 
  geom_boxplot()
pic2 = ggplot(coronal) + 
  aes(x = Region, fill = Fr) + 
  geom_bar()
pic3 = ggplot(int.df, aes(dB, Dauer)) + 
  geom_line() + 
  geom_point()
grid.arrange(pic1, pic2, pic3,  ncol=3, nrow =1)
```

# 11. Weitere Modifizierungen
```{r,eval=FALSE}
help(theme)
```

# 12. Zwei Unabhängige Variablen

In der Regel haben wir es mit einer abhängigen Variable zu tun,
die vermutlich von mindestens einer anderen Variable (die wir dann unabhängige Variable(n) nennen) abhängig ist (daher der Name "abhängige Variable"), d.h. sich durch deren Einfluss verändert (so vermuten wir zumindest).

Bislang haben wir Abbildungen erzeugt, die jeweils darstellen sollten, wie eine Variable sich durch den Einfluss einer anderen Variable ändert. Wir können aber auch abbilden, wie die abhängige Variable durch zwei oder drei unabhängige Variablen beeinflusst wird. Hierzu verteilen wir eine der unabhängigen Variablen wie bislang auf die x-Achse (entweder kontinuierlich, oder als stufiger, kategorialer Faktor).

```{r}
head(asp)
```

Nehmen wir nun an, wir wollen wissen, wie die Dauerwerte in `d` in Abhängigkeit von sowohl `Kons` als auch von `Bet` variiert. Die Dauerwerte sind kontinuierlich und numerisch, sowohl die Konsonantenkategorien als auch die Betonung (unbetont vs. betont) sind kategorial. Für eine numerische abhängige Variable und eine kategoriale unabhängige Variable benutzen wir eine Boxplot-Darstellung:

```{r}
ggplot(asp) + 
  aes(y = d, x = Kons) + 
  geom_boxplot() 

ggplot(asp) + 
  aes(y = d, x = Bet) + 
  geom_boxplot() 
```

Wie man sieht, beschreiben beide Abbildungen jeweils auf der y-Achse die Werte für `d`, wenn eine der jeweils zwei Bedingungen auf der x-Achse gegeben ist. Bezeichnen wir also die Abhängigkeit einer Variable y durch eine Variable x mit

`y ~ x`

Einen solchen Zusammenhang kann man sehr gut in einer zweidimensionalen Ebene abbilden, da wir es eben mit zwei Dimensionen zu tun haben. Wir haben es aber hier mit zwei unabhängigen Variablen und damit also mit drei Dimensionen zu tun: 

`y ~ x1 + x2`

Wir haben nun zwei Möglichkeiten: Wir können die zusätzliche Variable:

- i. mit Farbe kodieren
- ii. getrennt in einem anderen Panel darstellen.

### Farbe (`col = x2` oder `fill = x2` in `aes()`)

```{r}
# y ~ x1 + x2
ggplot(asp) + 
  aes(y = d, x = Kons, col = Bet) + 
  geom_boxplot() 
# kann auch y ~ x2 + x1 sein
ggplot(asp) + 
  aes(y = d, x = Bet, col = Kons) + 
  geom_boxplot() 
# statt der Umrissfarbe kann man die Füllfarbe mit fill setzen
#HINWEIS: keinesfalls col UND fill gleichzeitig benutzen
ggplot(asp) + 
  aes(y = d, x = Bet, fill = Kons) + 
  geom_boxplot() 
```

### Linientyp (`lty = x2` in `aes()`)

```{r}
# y ~ x1 + x2
ggplot(asp) + 
  aes(y = d, x = Kons, lty = Bet) + 
  geom_boxplot() 
```

### Panels (`facet_wrap(~x2)`)

```{r}
# Zu ii. Getrennte Panels, x1 in einem Panel, x2 in dem anderen mit der Funktion facet_wrap()
ggplot(asp) + 
  aes(y = d, x = Kons) + 
  geom_boxplot() + 
  facet_wrap(~Bet)
```

Was ist der Einfluss des Dialektes (Region) und Sozialklasse (Socialclass) auf die Artikulationsstelle (Fr) in dem Data-Frame `coronal`?

- y (Fr): kategorial
- x1 (Region): kategorial
- x2 (Sozialklasse): kategorial


Alle Variablen in `y ~ x1 + x2` sind kategorial: daher `geom_bar()`:

```{r}
# Nur x1 (also ohne zwischen den Sozialklassen zu differenzieren)
p1 = ggplot(coronal) + 
  aes(fill = Fr, x = Region) + 
  geom_bar() 
p1
# Mit x2
p2 = facet_wrap(~Socialclass)
p1 + p2

# Mit eigenen Farben
farben = c("cyan", "gold")
p3 = scale_fill_manual(values = farben) 
p1 + p2 + p3
```

Inwiefern wird die Dauer (Dauer) von der Intensität (dB) - getrennt pro Person (Vpn) - im Data-Frame int.df beeinflusst?

- y: Dauer (numerisch)
- x1: Intensität (numerisch)
- x2: Vpn (kategorial)

In `y ~ x1 + x2` sind zwei Variablen (eine davon die abhängige) numerisch. Daher benutzen wir `geom_point()`:

```{r}
ggplot(int.df) + 
  aes(y = dB, x = Dauer, col = Vpn) + 
  geom_point() 
```

Oder:

```{r}
# Getrennte Panels, x1 in einem Panel, x2 in dem anderen.
ggplot(int.df) + 
  aes(y = dB, x = Dauer) + 
  geom_point() + 
  facet_wrap(~Vpn)
```

Selbstverständlich geht auch eine Kombination von farblicher Kodierung und getrennten Panels:

```{r}
ggplot(int.df) + 
  aes(y = dB, x = Dauer, col = Vpn) + 
  geom_point() + 
  facet_wrap(~Vpn)
```


### Exkurs: Reihenfolge von Faktorstufen ändern

Bei `coronal` waren die Faktorstufen von `Socialclass` in einer bestimmten Reihenfolge, nämlich in alphabetischer:

```{r}
levels(coronal$Socialclass)
```

Hierbei steht aber `LM` für "lower middle class", `UM` für "upper middle class", und `W` für "working class". Es wäre stattdessen also die aufsteigende Reihenfolge `W`, `LM`, `UM` wünschenswert.

Um das zu ändern, können wir (in diesem Fall) `relevel()` verwenden:

```{r}
# in alphabetischer, aber eigentlich "falscher" Reihenfolge:
ggplot(coronal) + 
  aes(fill = Fr, x = Socialclass) + 
  geom_bar(position = "fill") +
  facet_wrap(~Region)

coronal$Socialclass = relevel(coronal$Socialclass,ref="W")
levels(coronal$Socialclass)
# und somit
ggplot(coronal) + 
  aes(fill = Fr, x = Socialclass) + 
  geom_bar(position = "fill") +
  facet_wrap(~Region)
```

Mit `relevel()` können Sie aber nur bestimmen, dass eine Faktorstufe, hier `W`, als neue Referenzklasse gewählt wird (daher `relevel(,ref = ...)`). Die anderen beiden Klassen bleiben in ihrer Reihenfolge unverändert, was hier zufällig mit unserer Wunschreihenfolge übereinstimmt. Dies wird aber sicher nicht immer der Fall sein.

Besser ist es daher, `factor()` mit dem Parameter `factor(,levels = ...)` zu benutzen, dann können wir wirklich die Reihenfolge ganz nach unseren Wünschen gestalten.

Stellen wir uns z.B. vor, dass wir die Reihenfolge diesmal absteigend (beginnend mit "upper middle class", endend mit "working class") haben wollen:

```{r}
levels(coronal$Socialclass)
###
coronal$Socialclass = factor(coronal$Socialclass, levels = c("UM","LM","W"))
#dann ist die Reihenfolge umgedreht
levels(coronal$Socialclass)
ggplot(coronal) + 
  aes(fill = Fr, x = Socialclass) + 
  geom_bar(position = "fill") +
  facet_wrap(~Region)
```