1 Organisatorisches

1.1 Software ↑

1.2 Kurs ↑

  • 6 Termine wöchentlich donnerstags 13.00-16.15:
    • 26.10. vor Ort in Biel
    • 2.11. Remote
    • 9.10. Remote
    • 16.11. vor Ort in Biel
    • 23.10. Remote
    • 30.11. vor Ort in Biel
  • ggf. Mittagspause kürzer und bis 15.45?
  • Corona Screencasts
  • Erreichbarkeit ausserhalb des Kurses
    • Forum
    • rudolf.farys(at)soz.unibe.ch
    • lukas.hobi(at)bfh.ch
    • olivier.lehmann(at)bfh.ch

1.3 Kompetenznachweis ↑

  • kleinere Übungen während der Kurstage (nicht bewertet aber sollen gelöst werden)
  • Note: Take-Home-Exam am Ende des Kurses (Einreichung bis 5. Januar 2024)

1.4 Kursinhalte ↑

  • Grundlagen, Workflow
  • Datenmanagement
  • Datenmanagement, Überblick über Funktionen zur Datenanalyse

2 Nützliches ↑

3 Warum R? ↑

4 Datenhandling und Objektorientierung ↑

5 RStudio

6 Zwei Fallbeispiele

6.1 Beispiel 1: Visualisierung eines Bewegungsprofils (Google Location API / JSON)

Rohdaten

{
  "locations" : [ {
    "timestampMs" : "1417727161925",
    "latitudeE7" : 469345574,
    "longitudeE7" : 74339415,
    "accuracy" : 31
  }, {
    "timestampMs" : "1417727100712",
    "latitudeE7" : 469345603,
    "longitudeE7" : 74339025,
    "accuracy" : 32
  }, {
    "timestampMs" : "1417727030427",
    "latitudeE7" : 469345829,
    "longitudeE7" : 74339350,
    "accuracy" : 30
  }, {
    "timestampMs" : "1417726978911",
    ...

Einlesen in R

library(jsonlite)
library(plyr)

raw <- fromJSON('LocationHistory.json')

Datenansicht in R

> head(raw$locations)
    timestampMs latitudeE7 longitudeE7 accuracy velocity altitude activitys heading
1 1446466295025  469344810    74340787       20       NA       NA      NULL      NA
2 1446466154216  469344810    74340787       20       NA       NA      NULL      NA
3 1446466096080  469344859    74340783       20       NA       NA      NULL      NA
4 1446465976033  469344859    74340783       20       NA       NA      NULL      NA
5 1446465853219  469344859    74340783       20       NA       NA      NULL      NA

> lapply(locs,class)
$timestampMs
[1] "character"

$latitudeE7
[1] "integer"
...

# -> kein ordentliches Datenobjekt sondern eine riesige Liste
# -> Datentypen/Skalen noch nicht richtig

Aufbereitung

locs = raw$locations

# Einen data.frame bauen, der die Daten enthält:

ldf = data.frame(t=rep(0,nrow(locs)))

# Zeitformat ist als Character (POSIX) in ms gespeichert. 
# Umwandeln in eine Zahl (in Sekunden)
ldf$t = as.numeric(locs$timestampMs)/1000
class(ldf$t) = 'POSIXct'

# lat/lon sind gemessen in "10^7""
ldf$lat = as.numeric(locs$latitudeE7/1E7)
ldf$lon = as.numeric(locs$longitudeE7/1E7)

# zum Spass noch ein paar Variablen bauen:
ldf$Jahr = substring(as.Date(ldf$t),1,4)
ldf$Sommer = substring(as.Date(ldf$t),6,7)%in%c("06","07","08","09")

fertige Daten

> head(ldf)
                    t      lat      lon Jahr Sommer
1 2015-11-02 13:11:35 46.93448 7.434079 2015  FALSE
2 2015-11-02 13:09:14 46.93448 7.434079 2015  FALSE
3 2015-11-02 13:08:16 46.93449 7.434078 2015  FALSE
4 2015-11-02 13:06:16 46.93449 7.434078 2015  FALSE
5 2015-11-02 13:04:13 46.93449 7.434078 2015  FALSE
6 2015-11-02 13:03:23 46.93449 7.434096 2015  FALSE

Plotten

bern = get_map('Bern, CH',13,scale=2,source="google",maptype="hybrid")
ldf.bern = mapclip(ldf,bern)
png("bernmap.png", width=1200, height=1200)
ggmap(bern) + geom_point(data=ldf.bern,aes(x=lon,y=lat,colour=Sommer),size=2,alpha=0.1)
dev.off()

6.2 Beispiel 2: Modell für Fussballwetten (football-data.co.uk)

  • Ziele
    • Ein Ratingsystem für Fussballmannschaften entwickeln
    • Ein Prognosemodell für Fussballwetten trainieren
  • Daten: Historische Ergebnisse und Wettquoten europäischer Mannschaften: http://www.football-data.co.uk/data.php
  • Schritte:
    • Kommaseparierte Daten automatisiert herunterladen (Crawler)
    • alle Daten poolen
    • aufbereiten
    • Rating-Algorithmus anwenden
    • Modelle bauen / Resultate ausgeben lassen
# für jede Website die csv Links grabben
for(url in urls) {
doc <- read_html(url)
links <- html_nodes(doc, "a") %>% html_attr("href")
csvlinks <- paste0("https://www.football-data.co.uk/", 
                   links[grep(pattern = ".csv", 
                              x = links)])

# Was bedeuten die Variablen?
# Div = League Division
# Date = Match Date (dd/mm/yy)
# HomeTeam = Home Team
# AwayTeam = Away Team
# FTHG = Full Time Home Team Goals
# FTAG = Full Time Away Team Goals
# FTR = Full Time Result (H=Home Win, D=Draw, A=Away Win)

# Alle Links parsen und zusammentackern
for(i in 1:length(csvlinks)) {
  tempdata <- as.data.table(read.csv(csvlinks[i]))
  try(dsocc <- rbind(dsocc,tempdata[, c("Div", "Date", "HomeTeam", "AwayTeam", "FTHG",
                                       "FTAG", "BbAvH", "BbAvD", "BbAvA"),
                                    with=FALSE], fill=TRUE, use.names=TRUE))
}
}

# Variablen definieren
dsocc[,win:=ifelse(FTHG>FTAG,1,ifelse(FTHG<FTAG,0,0.5))] # Sieg
dsocc[,numDate:=as.numeric(as.Date(Date,format="%d/%m/%y"))] # Datum numerich
dsocc[,HomeTeam:=as.character(HomeTeam)] # Datentyp fixen
dsocc[,AwayTeam:=as.character(AwayTeam)] # Datentyp fixen

Algorithmus (vereinfachter Ausschnitt)

  # Ähnlich wie die Elo-Zahl 
  # bekannt durch Schach, mittlerweile in Variationen aber sehr verbreitet:

  # 1. Startwerte setzen für alle Mannschaften
  # 2. Beim ersten Spiel im Datensatz anfangen -> Algorithmus anwenden, Ratings updaten
  # 3. Durch die Daten Zeile für Zeile gehen und Algorithmus anwenden, Ratings updaten
  # 4. Historie der Ratings (Mannschaftsstärken) speichern/mitprotokollieren

  # Kernstück:

  # Erwartungswerte für Tore und Gegentore
  E_a <- lambda * (off_a + h1) / def_b
  E_b <- lambda * (off_b - h2) / def_a
  
  # Tatsächlicher Spielstand (Tore/Gegentore):
  S_a <- dbl$FTHG[i]
  S_b <- dbl$FTAG[i]
  
  # neues Rating definieren auf Basis vom bisherigen Rating und den 
  # Abweichungen von der Erwartung
  neu_off_a <- off_a + K * (S_a - E_a)
  neu_def_a <- def_a - K * (S_b - E_b)
  neu_off_b <- off_b + K * (S_b - E_b)
  neu_def_b <- def_b - K * (S_a - E_a)
off def team mean
4162.644 3462.252 Bayern Munich 3812.448
3537.812 2576.271 Dortmund 3057.041
2946.873 2627.050 Leverkusen 2786.962
2459.675 2081.035 M’gladbach 2270.355
2453.434 2071.038 Hoffenheim 2262.236
2303.530 2368.270 Ein Frankfurt 2335.900
2294.210 2460.631 Wolfsburg 2377.421
2220.462 2165.347 Stuttgart 2192.904
2049.700 2073.904 Mainz 2061.802
2047.867 1914.318 Werder Bremen 1981.092
1949.648 2159.870 FC Koln 2054.759
1881.011 1916.348 Augsburg 1898.680
1804.919 1821.756 Hertha 1813.337
1756.083 2105.402 Ingolstadt 1930.742
1752.021 1788.077 Darmstadt 1770.049
1656.994 1790.680 Hannover 1723.837
1562.714 1740.381 Schalke 04 1651.548
1499.047 2073.383 Hamburg 1786.215

7 Arbeitsverzeichnis, Objekte und Workspace ↑

7.1 Beispiel: Arbeitsverzeichnis, Objekte und Workspace ↑

# Kommentare beginnen mit #
# Alles in der Zeile nach # wird von R ignoriert

# 1. R starten
# 2. Arbeitsverzeichnis (mit Schreibrecht) anlegen 
#    und R auf dieses Verzeichnis setzen: setwd()
# 3. Code ausführen mit STRG + R (Oder Icon "Run" oben rechts)

5 + 5

getwd() # Arbeitsverzeichnis anzeigen
# setwd() # Arbeitsverzeichnis definieren
setwd("meinpfad")

dir() # Arbeitsverzeichnis anzeigen

a <- 50 # Erzeugt Objekt a (Vektor der Länge 1) mit dem einzelnen Wert 50
a

# Objekt erzeugen, dass die Vornamen der Beatles enthält
die.beatles <- c("John", "Paul", "George", "Ringo")
die.beatles

# Mit c() - concatenate lässt sich auch ein Zahlenvektor bauen:
b <- c(1, 2, 3, 4)

# oder kürzer
b <- seq(1, 4)

# oder noch kürzer
b <- 1:4

# Objektnamen dürfen keine Leerzeichen haben. Ferner empfiehlt es sich - und _ zu meiden 
# siehe Google R Style Guide
# Namen sollten aussagekräftig sein. Namensgebung sollte im ganzen Code-File konsistent 
# sein (Punkte, Gross-/Kleinschreibung)

ls() # Workspace anzeigen

# Objekte a, b und die.beatles speichern in "beispiel1.RData"
save(a, b, die.beatles, file = "beispiel1.RData")

# RData Format ist sehr advanced (minimaler Speicherplatz und schnelle Schreib- 
# und Ladegeschwindigkeit) und empfiehlt sich immer, sofern der Export in andere 
# Formate (z.B. csv) nicht dringend nötig ist.

# Objekte löschen
rm(a, b, die.beatles) #oder:
rm(list=ls()) # löscht den gesamten Workspace
ls() # was ist jetzt noch im Workspace?
die.beatles # nicht mehr da

load("beispiel1.RData") # gespeichertes Objekt laden

die.beatles # wieder da

8 Packages ↑

8.1 Beispiel: Packages ↑

# nehmen wir an wir möchten eine SPSS Datei einlesen:

??spss

# liefert Hinweise auf die Funktionen read.spss und read_por aus den Paketen
# foreign und haven 
# allerdings heutzutage einfacher, zu googeln!

install.packages("foreign") # Install package

?read.spss 
# geht erst wenn das Package auch geladen ist

library(foreign)
?read.spss

9 Arithmetische- und logische Operatoren ↑

9.1 Beispiel: Rechnen und Vergleichen ↑

# Rechnen
ergebnis <- (23+24)*11/(18+15)*5
ergebnis

# Funktionen
log(2) 
cos(2)

# Vergleich
x <- -3:3
x

# sind die Elemente von x gleich 0?
x == 0

# grösser 0?
x > 0

# kleiner 0?
x < 0

# grösser gleich 0?
x >= 0

# kleiner gleich 0?
x <= 0

# ungleich 0?
x != 0

# grösser als -1 aber kleiner als 1
x > -1 & x < 1

# grösser als  1 und kleiner als -1
x > 1 & x < -1

# grösser als 1 oder kleiner als -1
x > 1 | x < -1

9.2 Übung: Rechnen und Vergleichen ↑

  1. Arbeiten Sie in einem Scriptfile innerhalb von RStudio. Bewahren Sie den Code für diese Übung (und generell für die weiteren Übungen) in diesem Scriptfile auf.
  2. Berechnen Sie folgende Terme in R
    1. \((3 + 4)^{2}\)
    2. \(\frac{-99}{33} + 42\)
    3. \(log(1)\)
    4. \((\sqrt{2})^{2}\)
  3. Prüfen Sie folgende Vergleiche:
    1. \(5 = 7\)
    2. \(5 \times 5 \geq 6 \times 4\)
    3. \(\sqrt{3} \neq cos(17)\)
  4. Für die mean() Funktion kann ein zusätzlicher Parameter trim angegeben werden. Finden Sie heraus, was mit diesem Parameter getan werden kann indem Sie die Hilfe von mean() lesen und halten Sie diese Info in ihrem Scriptfile fest.

9.3 Lösung: Rechnen und Vergleichen ↑

# 2.
(3 + 4)^2
-99/33 + 42
log(1)
sqrt(2)^2

# 3.
5 == 7
5 * 5 >= 6 * 4
sqrt(3) != cos(17)

# 4.
# help(mean)
# ?mean
# the fraction (0 to 0.5) of observations to be trimmed from each end of x 
# before the mean is computed.

10 Klassen/Datentypen ↑

10.1 Übersicht Datentypen ↑

# Homogen

# Integer Vektor
x <- 1:9
class(x)
x


# Numerischer Vektor
x <- c(1.3, 2.4, 3.5)
class(x)
typeof(x)
x


# Logischer Vektor
x <- -3:3
y <- x >= 0
class(y)
y

# String/Character Vektor
x <- c("a", "b", "c", "d", "f")
class(x)
x


# Matrix
matrix.2mal3 <- matrix(c(1,2,11,12,20,30), nrow = 2, ncol=3)
class(matrix.2mal3)
typeof(matrix.2mal3)
matrix.2mal3


# Array (z.B. 3 Matrizen)
arraybsp <- array(1:50, c(5,5,2)) # Zahlen von 1 bis 50 einem Array mit 2 5x5 Matrizen 
class(arraybsp)
arraybsp



# Heterogen

# Liste
liste <- list(a= c(4:8), b = c(1:3), c = c(2:10))
class(liste)
liste



# Data frame: Kann man anschauen durch Aufruf von swiss oder fix(swiss)
class(swiss)
swiss


# Faktoren
sex <- c(0,0,1,1)
factor(sex,labels=c("Mann","Frau"))

# warum ist das nuetzlich?
a <- rep(c("Haus", "Strasse"), 10^5)
object.size(a)
object.size(as.factor(a))
 
a<-rep(c("Haus", "Strasse"), 5)
object.size(a)
object.size(as.factor(a))

# Funktionen: z.B. cos(); mean()
class(mean)
mean

#  Datum
# https://stat.ethz.ch/R-manual/R-devel/library/base/html/strptime.html
daten.unformatiert<-c("19990123","20110324","20100412")
daten.als.daten<-as.Date(daten.unformatiert,format = "%Y%m%d")
tage.seit<-c(17000,14000,13000)
tage.als.daten<-as.Date(tage.seit,origin="1970-01-01")

10.2 Sonderwerte ↑

  • Inf und -Inf: Positiv und negativ unendlich
  • NaN: “Not a number”, z.B. 0/0
  • NA: fehlender Wert (Missing)
# Wichtiger Hinweis zu fehlenden Werten:
x <- c(1, 2, NA, 4)

#falsch:
x == NA
x == "NA"

#richtig:
is.na(x)

# siehe hierzu auch:
is.infinite(x)

10.3 Vektoroperationen ↑

1:100
1:100 * 3
1/100 * 3

# Vektor erzeugen und speichern
x <- 3 * 1:100
x

# oder über concatenate
x <- c(1, 2, 3, 4, 5)
x

# Vektor mit 15 Elementen. 5 Wiederholungen von 1,2,3
x <- rep(1:3, time=5) #default ist times: verkettet den Vektor 1 bis 3 5 mal
x
x <- rep(1:3, each=5) #each wiederholt jedes Element 5 mal
x
# Vektor mit Zahlen von 10 bis 100 in Zehnerschritten
x <- seq(10, 100, 10)
x

# Verknüpfen
x <- c(a = 10, b = 20, c = 30, d = 40)
x

10.4 Zugriff auf Datenbereiche ↑

  • Zugriff auf Elemente:
    • einzelne Elemente
      • durch Angabe ihrer Position im Vektor
      • oder durch einen Namen falls vorhanden
    • mehrere Elemente
      • Zugriff über logische Operatoren (TRUE/FALSE), d.h. über Bedingungen
      • oder eine Sequenz von Positionen oder Namen
x <- 11:200
x
x[1]              # erstes Element von x
x[1:10]           # die ersten 10 Elemente von x
x[-(11:100)]      # Alle Elements von x ausser die Positionen 11 bis 100

# Allgemeiner gibt es drei Möglichkeiten:
x <- c(a = 10, b = 20, c = 30, d= 40) # vergibt direkt Namen, alternativ:
# names(x) <- c("a", "b", "c", "d")

# Möglichkeit 1
x[1:2] # über die Position

# Möglichkeit 2
x[c("a", "c")] # über den Namen

# Möglichkeit 3
x[x < 20 | x >= 30] # über eine Bedingung


# head / tail
x <- 1:100
head(x)
tail(x)

10.4.1 Übung: Vektoren ↑

  1. Arbeiten Sie weiter in ihrem Übungs-Script-File
  2. Erzeugen Sie einen Vektor x der Länge 50, der die Zahlen von 1 bis 5 zehn mal wiederholt.
  3. Erzeugen Sie einen Vektor y (Länge 3), der die Elemente von x an den Positionen 12, 20 und 50 enthält.
  4. Erzeugen Sie einen Vektor freunde mit drei Namen Ihrer Wahl.

10.4.2 Lösung: Vektoren ↑

# 2.
x <- rep(1:5, 10)
x
length(x)

# 3.
y <- x[c(12, 20, 50)]
y

# 4. 
freunde <- c("Bernd", "Frank", "Franz")

10.5 Faktoren und Listen ↑

  • Faktoren
    • Class: "factor"
    • factor(): Erzeugt einen Faktor
    • sind ein Weg, kategoriale Daten (nominal/ordinal) sauber zu speichern
    • sind Vektoren mit Attributen
    • levels(): Zeigt die Kategorien eines Faktors
    • as.numeric(): Zwingt den Faktor in einen numerischen Vektor
  • Listen
    • Class: "list"
    • list(): Erzeugt eine Liste
    • Können beliebige Objekte/Objekttypen beinhalten
    • list$switzerland: Greift auf das Element switzerland der Liste list zu
    • list[2]: Zweites Element von list
    • Listen werden häufig als Rückgabewert von Funktionen verwendet, z.B. als Behältnis für Schätzergebnisse
  • Im Laufe des Kurses werden noch weitere Objektklassen auftauchen

10.5.1 Beispiel: Faktoren und Listen ↑

# Faktor erzeugen

sex <- factor(c(rep(0, 50), rep(1, 50)), labels = c("Mann", "Frau"))

# Hinweis zu as.numeric(): manchmal hat man folgendes Problem:
jahrgang <- factor(c("2000", "2000", "2001", "2002")) # falscher Datentyp "character"
# umwandeln in numerische Werte liefert aber nicht, was wir wollen
as.numeric(jahrgang)

# besser:
as.numeric(as.character(jahrgang))

# Eine Liste erzeugen:
kursteilnehmer <- list(Maenner = c("Simon", "Peter", "usw."),
                    Frauen = c("Daniela", "Johanna"))
kursteilnehmer

# Zuriff auf Elemente der Liste
kursteilnehmer$Maenner
kursteilnehmer[1]
kursteilnehmer[[1]]
kursteilnehmer["Maenner"]
kursteilnehmer[["Maenner"]]

# was ist der Unterschied?
class(kursteilnehmer[["Maenner"]])
class(kursteilnehmer["Maenner"])
# d.h. wenn man mit den Elementen aus einem Listenelement arbeiten will,
# braucht es doppelte Klammern oder Zugriff über "$"

length(kursteilnehmer)
kursteilnehmer$Frauen[2]
kursteilnehmer[[2]][2]

10.6 Data Frames ↑

  • Data frames sind das typische Format für Datensätze
  • …es handelt sich um eine Liste von Vektoren derselben Länge siehe Wickham
  • …ähnelt Matrizen aber die Spalten können unterschiedliche Datentypen beinhalten
  • data.frame(): erzeugt einen data frame
  • as.data.frame(): konvertiert in einen data frame
  • order(): sortiert Daten
  • summary() und str(): Überblick über data frames
  • head() and tail(): erste/letzte Zeilen inspizieren
  • names(): zeigt Spaltennamen
  • object$var1: Greift direkt auf die Spalte var1 im data frame object zu
  • na.omit(): Zeilenweise Ausschluss von fehlenden Werten, d.h. Zeilen die mindestens 1 Missing beinhalten

10.6.1 Beispiel: Data Frames ↑

# fertige Daten sind oft data frames:
titanic <- read.dta("http://www.stata-press.com/data/kkd/titanic2.dta")
is.data.frame(titanic)

# man kann sich auch leicht selber einen bauen
obst <- c("Apfel", "Apfel", "Birne")
gemuese <- c("Tomate", "Karotte", "Karotte")
id <- 1:3
df <- data.frame(id, obst, gemuese)
df

# Ansteuern von Zeilen und Spaltenpositionen
df$obst
df[, "obst"]
df[3, "gemuese"]
df[3, 3]

# d.h. es gibt wie bei Vektoren diverse Möglichkeiten: Position, Name ($ oder "") oder Bedingung

10.6.2 Spalten hinzufügen und löschen ↑

# To add a column to a data.frame you can also use the $
df$drinks <- c("Milch", "Cola", "Bier")
df

# To delete a column, assign a NULL:
df$gemuese <- NULL
df

10.6.3 Übung: Data Frames ↑

  1. Verwenden Sie den R-eigenen Datensatz swiss. Verschaffen Sie sich ggf. einen Überblick mit ?swiss.
  2. Lassen Sie sich die zehnte bis zwölfte Zeile des Datensatz ausgeben.
  3. Lassen Sie sich die Spalten Education und Catholic ausgeben, jedoch nur für Fälle, deren Kindersterblichkeit zwischen 20 und 22 liegt.

10.6.4 Lösung: Data Frames ↑

swiss[10:12, ]

swiss[swiss$Infant.Mortality > 20 & swiss$Infant.Mortality < 22, c("Education", "Catholic")]

11 Kontrollstrukturen und Schleifen ↑

11.1 Beispiel: Kontrollstrukturen und Schleifen ↑

# For Schleife
for (x in 1:10) {
  print(sqrt(x)) 
}

# aber besser:

sqrt(1:10) # da Funktionen in R i.d.R. sowieso vektorisiert arbeiten

# Schleifen machen aber Sinn, wenn es Abhängigkeiten zwischen den Durchläufen gibt:
x <- 0
for(i in 1:10) {
  x <- x+i
  print(x)
  }
# x wird immer weiter inkrementiert
# mit ein bisschen Überlegen gibt es aber auch hier eine Lösung die performant ist und vektorisiert arbeitet:
cumsum(1:10)

# es braucht also schon etwas kompliziertere Abhängigkeiten


# Sequenzen in For-Schleifen können auch Character sein:
namen <- c("Alfred", "Jakob", "Peter")
for (name in namen) {
  print(paste("Hallo", name))
}

# einfacher aber:
paste("Hallo", namen) # weil vektorisiert

#sinnvolleres Beispiel:

for (dataset in c("data1.csv", "data2.csv", "data3.csv")) {
  read.csv(pfad/dataset)
  # Anweisungen, 
  # z.B. Datenbereinigung, Appending (rbind), Modellschätzungen, etc.
}

# das Beispiel könnte aber auch mit Hilfer einer selbstgeschriebenen Funktion gut gelöst werden (siehe nächstes Kapitel)


# Durch Spalten loopen
for (column in 2:6) { # this loop runs through 2 to 6
        print(names(swiss)[column])
        print(mean(swiss[, column]))
}

# aber wieder geht es einfacher und schneller:

colMeans(swiss[, 2:6]) # oder
apply(swiss[, 2:6], 2, mean) 

# die 2 verweist auf "spaltenweise" (1 wäre zeilenweise). 
# D.h. für jede Spalte der Daten wird mean() angewendet



# While Schleife
x <- 0 # Startbedingung sollte gelten
while(x < 13) {
  x <- x+1 # inkrementieren, da sonst die Bedingung für immer gilt -> Endlosschleife
  print(x)
} 

# wird wiederholt solange x<13==TRUE
# sicherstellen, dass irgendwann das Kriterium FALSE wird!

# Beispiel für eine sinnvole while-Schleife: Abfrage einer Web-Ressource, die nicht immer erreichbar ist. while(keinen erfolg) try(ressource abfragen)


# if-Beispiel:

# Daten einlesen aus einer Liste von Files
setwd("C:/path/to/some/excel/files")
myfiles <- list.files()
# manche sind aber nun xls, und andere xlsx:

library(tools)
for(file in myfiles) {
  if(file_ext == "xls") {
    daten <- read.xls(file)
  }
  if(file_ext == "xlsx") {
    daten <- read.xlsx(file)
  }
}

# if prüft immer nur genau eine Bedinung. 


# es unterscheidet sich dadurch vom Kommando ifelse, das vektorisiert arbeitet
# ifelse Beispiel:

a<- sample(1:100, 10)
b<-ifelse(a < 50, "Nicht bestanden", "Bestanden")
b

# ifelse prüft einen Vektor von Bedingungen. Naheliegenderweise ist so
# ein Konstrukt also auch gut zur Datenaufbereitung geeignet.
# obiges Beispiel ist identisch mit aber einfacher als:

b[a < 50] <- "Nicht bestanden"
b[a >= 50] <- "Bestanden"
b
  • Achtung: Die meisten Schleifen sind vermeidbar und es gibt effizientere Lösungen, z.B. Funktionen aus dplyr, sapply()/lapply()/apply() oder spezielle Funktionen, z.B. rowMeans(), colMeans(), colSums().
  • Überlegen Sie, ob eine Schleife wirklich nötig ist für ihr Problem (sehr komplizierte und interdependente Probleme z.B. können oft nicht vektorisiert werden und Schleifen sind dann ggf. unumgänglich).
  • Fall Sie Schleifen verwenden: Ziehen Sie soviel Code wie möglich aus der Schleife raus (alles statische) um nicht unnötig Performance einzubüssen. Vermeiden Sie innerhalb der Schleifen komplexe Datentypen (data.frames) und Funktionen die damit zusammenhängen, z.B. cbind()/rbind(). Verwenden Sie besser Vektoren, Matrizen oder data.tables (kommt in den folgenden Kapiteln).

11.2 Übung: Kontrollstrukturen und Schleifen ↑

  1. Schreiben Sie eine Schleife die 10 Durchläufe hat. In jedem Durchlauf sollen 100 standardnormalverteilte Zufallszahlen gezogen werden (rnorm(100)). Berechnen Sie innerhalb der Schleife Mittelwert und Standardabweichung dieser 100 Werte (mean()/sd()) und geben Sie diese aus (print).
  2. Verwenden Sie für dieselbe Aufgabe nun die Funktion tapply. Erzeugen Sie hierfür zwei Vektoren. Der erste soll alle 1000 Zufallszahlen beinhalten, der zweite die Samplezugehörigkeit (Vektor der Länge 1000 mit 100 1ern, 100 2ern, …). Diese zwei Vektoren können dann an tapply “gefüttert” werden.

11.3 Lösung: Kontrollstrukturen und Schleifen ↑

# 1. 
for (i in 1:10) {
  x <- rnorm(100)
  mittel <- mean(x)
  std <- sd(x)
  print(paste("Mittelwert:", mittel, "Standardabweichung:", std))
}

# Alternative mit tapply:
x <- rnorm(1000)
sample <- rep(1:10, each=100)
tapply(x, sample, mean)
tapply(x, sample, sd)

12 Funktionen ↑

12.1 Beispiel: Funktionen ↑

print("Hallo")
print

lm

# Beispiel: Definieren einer Funktion
wurzel <- function(x) {
  x^0.5
}


wurzel(4)

# Namen eingeben um Inhalt zu sehen
wurzel

12.2 Übung: Funktionen ↑

  1. Schreiben Sie eine Funktion, die Input gemäss folgender Formel verarbeitet: \((x + 2)^{2}\) und wenden Sie diese Funktion auf einen Vektor von Zahlen (1 bis 10) an.
  2. Schreiben Sie eine Funktion mit zwei Inputvariablen x und w, die einen gewichteten Mittelwert berechnet (x gewichtet mit w) und zurückgibt: \(\frac{\sum_{i=1}^{N}w_{i}x_{i}}{\sum_{i=1}^{N}w_{i}}\). Testen Sie ihre Funktion mit x <- 1:5 und w <- c(2,4,5,6,7).
  3. Bonus: Gehen Sie zurück zur tapply-Aufgabe. Wie könnte eine Lösung auf Basis einer selbstgeschriebenen Funktion aussehen (konzeptionell, lauffähiger Code nicht nötig)?

12.3 Lösung: Funktionen ↑

# 1.
plus2hoch2 <- function(x) {
  (x+2)^2
}

plus2hoch2(1:10)

# 2. 
w <- c(2, 4, 5, 6, 7)
x <- 1:5
gew.mittelwert <- function(x, w) sum(w * x) / sum(w)
gew.mittelwert(x, w)

# 3.
drawSample <- function(i) {
  x <- rnorm(100)
  data.frame(mean.val = mean(x), std = sd(x))
}

do.call(rbind, lapply(1:10, drawSample))

# noch kompakter mit Package data.table:
library(data.table)
rbindlist(lapply(1:10, drawSample))

13 Code lesbarer schreiben mit “Piping” ↑

13.1 Beispiel: Code lesbarer schreiben mit “Piping” ↑

daten <- rnorm(100)

# normale Schreibweise
mean(daten)

# mit Piping
daten %>% mean(.)
# der Punkt sagt "hier kommen die Daten rein, die links von der Pipie stehen"

# kürzer
daten %>% mean()

# oder sogar
daten %>% mean

# derart lassen sich aber mehrere Schritte verketten, z.B.:

daten <- c("1", "2", "3") # Strings statt Zahlen

daten %>%
  as.numeric %>% # umwandeln
  plus2hoch2

# später mehr zum Thema Piping!

14 Kontingenztabellen und einfache Tests ↑

14.1 Beispiel: Kontingenztabellen und einfache Tests ↑

# zwei fiktive Vektoren erstellen
auto <- factor(c(1, 0, 0, 0, 0, 
          1, 0, 1, 0, 0, 
          1, 1, 0, 0, 1, 
          1, 1, 0, 0, 0), 
          labels = c("kein Auto", "Auto"))
geschlecht <- c(rep("Frau", 10), 
                rep("Mann", 10))

# Tabelle
table(geschlecht)
table(geschlecht, auto)

# Chi2-Test
chisq.test(table(geschlecht, auto))

# Tabelle in Prozent
100*prop.table(table(geschlecht))
100*prop.table(table(geschlecht, auto))
round(100*prop.table(table(geschlecht, auto)), 2) # gerundete Werte
round(100*prop.table(table(geschlecht, auto), margin=1), 2) # margin=1 berechnet Zeilenprozente, margin=2 Spaltenprozente

# oder mit Piping
table(geschlecht, auto) %>%
  prop.table(., margin = 1) %>%
  round(., 2) * 100

# Lagemasse

## Mean
mean(x)

## Median
median(x) 
sort(x)

## gibt es verkürzt über die generische Funktion summary
summary(x)

# Streuung, z.B.

sd(x)
var(x)

# Korrelation zwischen Vektoren
cor(x, y)
cor(x, y, method = "spearman") # Rangkorrelation
cov(x, y)

# Mittelwertvergleich
t.test(x, y)

14.2 Übung: Kontingenztabellen und einfache Tests ↑

  1. Laden (ggf. vorher installieren) Sie das Paket foreign. Laden Sie die Daten über den Untergang der Titanic mit folgendem Befehl: titanic <- read.dta("http://www.stata-press.com/data/kkd/titanic2.dta") und machen Sie sich ein bisschen mit den Daten vertraut (z.B. head() oder summary()).
  2. Erzeugen Sie eine Kreuztabelle mit den Variablen class und survived.
  3. Berechnen Sie Median und Mittelwert des Alters für Überlebte und Gestorbene.
  4. Bonus: Unterscheiden sich die beiden Gruppen bzgl. ihres Alters? Was sagt uns hier ein t-Test?

14.3 Lösung: Kontingenztabellen und einfache Tests ↑

### 1.
library(foreign)
titanic <- read.dta("http://www.stata-press.com/data/kkd/titanic2.dta")
head(titanic)
summary(titanic)

### 2. 
# survived aus Konvention auf die Y-Achse, class auf die X-Achse (Y ist die abhängige, X die erklärende Variable)
table(titanic$survived, titanic$class) # Faktor crew ist hier nicht sauber gelabelt

### 3.
mean(titanic$age2[titanic$survived == "yes"])
mean(titanic$age2[titanic$survived == "no"])

# alternativ mit tapply
tapply(titanic$age2, titanic$survived, mean)

# ja sie unterscheiden sich um ca. 5 Jahre. Die Verstorbenen sind älter. Man könnte an der Stelle auch noch einen Signifikanztest machen:
t.test(titanic$age2 ~ titanic$survived)
# dieser zeigt, dass die 5 Jahre nicht zufällig sondern systematisch sind.

15 Workflow ↑

16 Mit Strings arbeiten

16.1 Beispiel: Mit Strings arbeiten ↑

# Gross-/Kleinschreibung
"trump" == "trump"
"Trump" == "trump"
tolower("Trump") == tolower("trump") 


# Trimmen
s <- "  Hello!   "
trimws(s)
trimws(s, "right")
trimws(s, "left")


# Extrahieren
s <- "I will build a great, great wall on our southern border, and I will have Mexico pay for that wall. Mark my words."
substr(s, 3, 6)
substr(s, 74, 79)


# Splitten
s <- c("To be blunt, people would vote for me. They just would.")
strsplit(s, ".", fixed = TRUE)


# Joinen
paste("one", "two", "three")
paste("one", "two", "three", sep = "-")

# paste erwartet normalerweise eine Reihe von ein-elementigen Inputs, nicht Vektoren:
s <- c("one", "two", "three")
paste(s)

# wir können paste() aber sagen es soll diesen Vektor "collapsen"
s <- c("one", "two", "three")
paste(s, collapse = " + ")
paste(s, collapse = "")


# Suchen und ersetzen
sub('Trump', 'Donald', 'Trump became president. Trump makes everyone great again.')
gsub('Trump', 'Donald', 'Trump became president. Trump makes everyone great again.')

16.2 Übung: Mit Strings arbeiten ↑

  1. Erzeugen Sie einen Vektor mit 100 Strings mit den Werten “objekt_1”, “objekt_2”,…

Verwenden Sie nun folgenden String:

s <- c("Kontostand: 100 EUR", "Kontostand: 150 EUR", "Kontostand: 185 EUR")
  1. Ersetzen Sie EUR durch CHF.
  2. Splitten Sie den String um nur den Geldbetrag zu erhalten

16.3 Lösung: Mit Strings arbeiten ↑

# 1. 
paste("objekt", 1:100, sep = "_")

# 2. 
s <- c("Kontostand: 100 EUR", "Kontostand: 150 EUR", "Kontostand: 185 EUR")
sub("EUR", "CHF", s)

# 3. 
strsplit(s, ": ", fixed = TRUE)
# um die gesplitteten Vektoren einzeln wieder einzusammeln kann man bspw. sapply() verwenden:
liste <- strsplit(s, ": ", fixed = TRUE)
sapply(liste, "[", 2)

# Man könnte so auch direkt an den numerischen Wert kommen, z.B.
liste <- strsplit(s, " ", fixed = TRUE)
as.numeric(sapply(liste, "[", 2))

# dafür gibt es aber effiztientere Lösungen, z.B. mit regulären Ausdrücken

16.4 Reguläre Ausdrücke / Pattern Matching

  • Einige R-Funktionen erlauben die Verwendung von sogenannten regulären Ausdrucken, beispielsweise sub/gsub

  • Reguläre Ausdrücke erlauben “wildcards” und andere Elemente um ein Suchmuster zu definieren

  • Eine Auswahl von Metacharacters:

    • Wildcards
      • .: matcht einen beliebigen Character.
      • [a-z]: lowercase a bis z
      • [A-Z]: uppercase A bis Z
      • [a-zA-Z]: lower a bis z und uppercase A bis Z
      • [0-9]: Ziffern
    • Quantifizierer
      • *: matcht mindestens 0 mal.
      • +: matcht mindestens 1 mal.
      • ?: matcht höchstens 1 mal.
      • {n}: matcht genau n mal.
      • {n,}: matcht mindestens n mal.
      • {n,m}: matcht zwischen n und m mal.
    • Position
      • ^: matcht am Anfang des Strings.
      • $: matcht am Ende des Strings.
    • Auslesen bestimmter Teile
      • ()

  • Es gibt viele weitere Metacharacter, siehe z.B.: https://stat.ethz.ch/R-manual/R-devel/library/base/html/regex.html

  • Siehe auch https://emailregex.com/ für ein ausführliches Regex Beispiel

16.4.1 Beispiel: Reguläre Ausdrücke / Pattern Matching ↑

# wildcards

s <- c("Trump", "Trrump", "Trrrump", "Tump")

grepl("T.ump", s) 
grepl("Tr+ump", s) 
grepl("Tr{1,2}ump", s) 
grepl("Tr{2,}ump", s) 
grepl("Tr*ump", s) 

grepl("T.*ump", s)

grepl("T[a-z]ump", c("Tump", "Trump", "Trrump", "Tdump"))

# ersetzen
gsub("T.*ump", "Donald", c("Tump", "Trump", "Trrump", "Tdump"))

# wir können Teile des Strings als Gruppe erfassen um diese zu verwenden

# z.B.:
gsub("(Donald) (Trump) (.*)", "\\2 \\1 \\3", "Donald Trump is president.")
# gsub erinnert sich hier an die drei Gruppen und gibt sie in der Reihenfolge 2., 1., 3. zurück

# Das Package stringr hat zudem einige convenience Funktionen die auf sub/gsub/grep/grepl aufbauen
library(stringr)
library(magrittr)

# CHF Beträge aus Text extrahieren ({1,5} -> min. 1 Ziffer, max. 5 Ziffern)
str_extract(c("Ein Brot kostet 4 CHF.", "Ein Auto kostet 50000 CHF."), "[0-9]{1,5} CHF")

# Strings matchen aber nur Teile davon zurückgeben:
str_match(c("Ein Brot kostet 4 CHF.", "Ein Auto kostet 50000 CHF."), "([0-9]{1,5}) CHF")
# die () definieren den Teil der zusätzlich zurückgegeben wird.

# Beispiel wie man an die zweite Spalte kommt:
str_match(c("Ein Brot kostet 4 CHF.", "Ein Auto kostet 50000 CHF."), "([0-9]{1,5}) CHF")[, 2] %>% as.numeric

16.4.2 Übung: Reguläre Ausdrücke / Pattern Matching ↑

Aus folgendem String

s <- "George Bush, born July 6, 1946, is a former American president. His full name is George Walker Bush."
  1. Extrahieren Sie George Bush’s Geburtsjahr und
  2. Seinen zweiten Vornamen

16.4.3 Lösung: Reguläre Ausdrücke / Pattern Matching ↑

# andere Lösungen ebenso möglich!
# 1.
gsub(".*([0-9]{4}).*", "\\1", s)

# oder mit stringr
library(stringr)
str_extract(s, "[0-9]{4}")

# 2.
gsub(".*George ([a-zA-Z]+) Bush.*", "\\1", s)

# oder mit str_match
str_match(s, "George ([a-zA-Z]+) Bush") # Buchstaben zwischen George und Bush
str_match(s, "George ([A-Z]{1}[a-z]+) Bush") # noch etwas expliziter mit grossen Anfangsbuchstaben

17 Datenmanagement ↑

17.1 Package dplyr ↑

  • Eines der populärsten Packages für data wrangling: dplyr von Hadley Wickham/Romain Francois bietet ein Toolset zur Datenaufbereitung
  • Siehe die dplyr Vignette und das Data Wrangling Cheat Sheet für einen sehr guten Überblick
    • filter(): wählt ein Subset von Zeilen (siehe auch slice())
    • arrange(): sortiert
    • select(): wählt Spalten
    • mutate(): erzeugt neue Spalten
    • summarise(): aggregiert (collapses) Daten zu einzelnen Datenpunkten
    • group_by(): Definiert Untergruppen in den Daten, damit o.g. Funktionen separat pro Gruppe angewandt werden können
    • uvvm.
    • dplyr kann sehr gut zusammen mit Piping verwendet werden, d.h. das Datenobjekt wird von Funktion zu Funktion weitergereicht durch %>% wodurch der Code deutlich besser lesbar und kompakter wird. (Base-R Pipe |>)

17.1.1 Beispiel: Package dplyr ↑

titanic <- read.dta("http://www.stata-press.com/data/kkd/titanic2.dta")

# install.packages("dplyr")
library(dplyr)

filter(titanic, class == "1st class", age2 < 18)

# konventionell wäre das komplizierter:
titanic[titanic$class == "1st class" & titanic$age2 < 18, ]

# zusätzlich Spalten selektieren:
titanic %>%
  filter(class == "1st class", age2 < 18) %>%
  select(sex, age2, survived)

# neue Variable "child" bauen
titanic %>%
  mutate(child = age2 < 18) %>%
  head()

# Auszählen wer gestorben ist nach Geschlecht und Kind (ja/nein)
titanic %>%
  mutate(child = ifelse(age2 < 18, "yes", "no")) %>%
  group_by(sex, child, survived) %>%
  summarise(n=n()) %>%
  arrange(sex, child, survived)

17.1.2 Übung: Package dplyr ↑

  1. Laden Sie die Daten http://www.farys.org/daten/allbus2008.dta (hierfür brauchen Sie die Funktion read.dta() aus dem Package foreign. Es handelt sich um Daten aus einer Bevölkerungsumfrage. R Version 3.4.x hat offenbar ein bislang nicht gefixtes Problem im Zusammenspiel mit read.dta(). Bei einem “factor level duplicated” Error bitte die Option convert.factors=FALSE verwenden.
  2. Wir brauchen nur die Variablen v151 (Geschlecht), v154 (Alter) und v386 (Einkommen). Legen Sie einen verkleinerten Datensatz an, der nur diese drei Variablen beinhaltet (und ordentliche Spaltennamen hat). Sie können in select direkt neue Spaltennamen vergeben: select(neuername = altername).
  3. Die Daten enthalten falsche Werte, die eigentlich als fehlende Werte deklariert sein müssten (z.B. 99997, 99998 und 99999 beim Einkommen und 999 beim Alter). Filtern Sie diese Werte weg (filter()).
  4. Berechnen Sie mit summarise() das Durchschnittseinkommen für jede Kombination aus Alter und Geschlecht (group_by(), summarise(), mean()).

17.1.3 Lösung: Package dplyr ↑

# Daten laden
allbus <- read.dta("http://www.farys.org/daten/allbus2008.dta")

allbus.agg <- allbus %>%
  select(geschlecht = v151, # die drei Variablen wählen und direkt umbenennen
         alter = v154,
         einkommen = v386) %>%
  filter(einkommen < 99997, alter < 999) %>% # fehlende Werte droppen
  group_by(geschlecht, alter) %>% # gruppieren
  summarise(m_einkommen = mean(einkommen)) # aggregieren

# Zusatz: Das ganze könnte man jetzt grafisch anschauen (müsste man ggf. etwas gröber gruppieren)
library(ggplot2)  
ggplot(allbus.agg, aes(x=alter,y=m_einkommen,color=geschlecht)) +
  geom_line()

17.2 Das Package data.table

  • Versucht dasselbe wie dplyr, ist aber stärker auf Performance und Funktionalität ausgerichtet während dplyr mehr die Lesbarkeit des Codes für den Menschen im Fokus hat
  • data.table bringt einen eigenen Objekttyp data.table mit, der mit allen data.frame Operationen abwärtskompatibel ist aber wesentliche Neuerungen beinhaltet:
  • Performance (schnellstes Paket für Datenoperationen, z.B. Sortieren, Matchen, Mergen, Gruppieren, Reshapen (dcast und melt ohne dass reshape2 geladen werden muss), Umbenennen von Spalten, Doubletten droppen, u.a.)
  • Memory efficient (beim Umbenennen/Umsortieren von Spalten oder Zeilen werden intern keine Kopien der Daten angelegt, die Daten sind lediglich “changed by reference”)
  • Generelle Form: DT[i,j,by]: “Nimm das DT an Zeilen i und berechne j, gruppiert nach by”
  • Kennt einige Konzepte zur eleganten Datenaufbereitung, etwa
    • .N: Anzahl Zeilen im Subset,
    • .BY: Platzhalter für den aktuellen Wert der Gruppierungsvariable vom Typ list. z.b. Nützlich, wenn man bei nach numerischen Werten gruppierten Variablenberechnungen auf den Gruppenwert zurückgreifen will. Oder wenn man bei nach kategorisch gruppierten Berechnungen auf den Namen zugreifen will (z.b. um zu ploten),
    • .SD: Enthält die Daten des Subsets das gerade verarbeitet wird (alle Spalten ausser by). Mit .SDcols spezifizieren sie, welche Spalten zu .SD gehören sollen. Nützlich bei der Anwendung einer Funktion auf mehrere, aber nicht alle Spalten eines DT.
    • .I: Index eines Falls innerhalb des Subsets,
    • .GRP: Group-Counter
  • Sehr sehr mächtig, aber Syntax ist etwas schwerer zugänglich und braucht mehr Eingewöhnung.
  • Beim Arbeiten mit grossen Datensätzen unerlässlich!

17.2.1 Beispiel: Package data.table ↑

library(foreign)
titanic <- read.dta("http://www.stata-press.com/data/kkd/titanic2.dta")
 
library(data.table)
titanic <- as.data.table(titanic) # als data.table definieren
titanic  # man beachte die andere Darstellung. data.table zeigt
# sicherheitshalber niemals alle Zeilen auf der Konsole
 
#Filtern ####
titanic[class=="1st class" & age2 < 18]
 
#Selektieren ####
titanic[,list(sex,class)]#shortcut für list(sex,class): .(sex,class)
titanic[,c("sex","class")]#data.table versteht c("name") als Liste
 
#Neue Variable berechnen ####
titanic[,child := ifelse(age2 < 18, "yes", "no")]# := um Variable by reference hinzuzufügen
titanic[,":="(child = ifelse(age2 < 18, "yes", "no"),
              oldperson = ifelse(age2 > 65, "yes", "no"))]#mehrere Variablen gleichzeitig
 
#Gruppierte Berechnungen ####
titanic[,Durchschnittsalter.Klasse:=mean(age2),by=.(class)] #Neue Spalte
titanic[,.(Durchschnittsalter.Klasse=mean(age2)),by=.(class)] #Aggregation
 
#Anwendung (eigener) Funktion auf eine Auswahl von Spalten####
# Z.b. Modusfunktion
Mode <- function(x) {
  val<-max(table(x))
  names(table(x))[table(x)==val]
}
 
#.SDcols wird als viertes Argument spezifiziert: definiert die zu bearbeitenden Spalten
vars <- c("survived","sex")
titanic[,lapply(.SD, Mode),,.SDcols=vars]
 
#Wenn man das Resultat in einer Variable behalten möchte: '(namensvektor):=' davor
varsnew <- paste0(vars,".Modus")
titanic[,(varsnew):=lapply(.SD, Mode),,.SDcols=vars]
 
#Chaining#####
titanic[, child := ifelse(age2 < 18, "yes", "no")][ # chaining: Äquivalent zum piping
  ,child.died := ifelse(child == "yes" &
                         survived == "no",
                       "yes", "no")]
 
#Tipps und Kniffs######
 
# Vorsicht beim Selektieren######
a<-c("sex","class")
titanic[,a] # Charaktervektoren können nicht als Objekte übergeben werden:
# data.table sucht immer nach Namen innerhalb des data.table, wenn es sie nicht findet
 
titanic[,..a] # Wenn Sie '..' vor einen Charactervector setzen funktioniert es.
titanic[,a,with=F] # Besser with=F (nicht within evaluieren).
 
# Dies ist nützlich, wenn sie nach Variablen mit bestimmtem Mustern suchen
# (z.b. mit grepl("s",names(titanic))). So versteht data.table auch logische
# Vektoren (oder Spaltennummern) als Selektionskriterien
 
titanic[,grepl("s",names(titanic)),with=F]
 
# Vorsicht beim Kopieren von Objekten#####
#Wenn man einen data.table kopiert, copy() verwenden, sonst werden neue Variablen auf das ursprüngliche Objekt geschrieben
 
titanic.copy <- copy(titanic)
 
#Beispiel für .I#####
a <- data.table(a = c(0,2,3,1,6,7,8))
a[, b := a[.I+1]]#Alternative zu lead() in dplyr
a[, c := c(NA, a[.I-1])]#Alternative zu lag() in dplyr
 
#neu auch
a[, b := shift(a, n = 1, type = "lag")]
a[, c := shift(a, n = 1, type = "lead")]

17.2.2 Übung: Package data.table ↑

  1. Laden Sie die Daten des European Social Survey (Auszug) mit load(url("http://www.farys.org/daten/ESS.RDATA")).
  2. Welches Land hat am meisten Befragte (Tipp: verwenden Sie .N und chainen sie die Bedingung [N==max(N)] um die Anzahl und den Namen des Landes auszugeben))?
  3. Welches ist das durchschnittlich glücklichste Land? Welches das unglücklichste?
  4. In welchem Land gibt es anteilsmässig am meisten “komplett glückliche Menschen” (d.h. der Anteil der Personen mit angegebener Zufriedenheitswert = 10 ist maximal; Tipp: berechnen Sie den mittleren Anteil pro Land und werfen sie den höchsten Wert aus)?
  5. Welches Land hatte von 2008 bis 2010 den grössten Rückgang im Durchschnittsglück (Tipp: über die eckigen Klammern können sie innerhalb des bearbeiteten Subsets Fälle auswählen, z.b. mean(happy[year==2008],na.rm=T))?
  6. Berechnen sie gleichzeitig den Modus bezüglich Glücklichkeit und Support für Umverteilung pro Land (Variable gincdiff 1=agree strongly,5=disagree strongly). Tipp: verwenden Sie dabei die Modusfunktion aus dem Beispiel oben.

17.2.3 Lösung: Package data.table ↑

library(data.table)
library(dplyr)
library(ggplot2)
load(url("http://www.farys.org/daten/ESS.RDATA"))
#gincdif: The government should take measures to reduce differences in income levels. 1=agree strongly,5=disagree strongly
#happy: Taking all things together, how happy would you say you are? #10 extremely happy 1 =Extremely unhappy
#uempla: Unemployed, actively seeking

#um Gruppenzähler kennenzulernen
#__________________________________

#2. Welches Land hat am meisten Befragte?#######
ess[,.N,Land][N==max(N)]

#dplyr
ess%>%
  group_by(Land)%>%
  count()%>%
  arrange(-n)

#Base
a<-table(ess$Land)
names(a)[a==max(a)] 
#Um Aggregatsfunktion kennenzulernen
#__________________________________

#3. Welches ist das durchschnittlich glücklichste Land? Welches das unglücklichste?#######
ess[,.(m.happy=mean(happy,na.rm=T)),Land][m.happy==max(m.happy)|m.happy==min(m.happy)]

#Alternativ: ohne order
hap <- ess[,.(m.happy=mean(happy,na.rm=TRUE)),Land]
setorder(hap,m.happy)
hap

#dplyr

ess%>%
  group_by(Land)%>%
  summarise(m.happy=mean(happy,na.rm=T))%>%
  arrange(-m.happy)


#Base R
ess<-ess[!is.na(ess$happy)&!is.na(ess$Land),]
agghapp <- aggregate(ess$happy,list(ess$Land),FUN=mean)
agghapp[order(agghapp$x),] 

#um Variablengenerator kennenzulernen (:=)
#__________________________________

#4. In welchem Land gibt es am meisten komplett glückliche Menschen?#######
ess[,":="(komplett=ifelse(happy==10,1,0))][
  ,.(`Anteil komplett glücklich`=mean(komplett,na.rm = T)),by=Land][order(`Anteil komplett glücklich`)]


#dplyr
ess%>%
  mutate(komplett=happy==10)%>%
  group_by(Land)%>%
  summarise(`Anteil komplett glücklich`=mean(komplett,na.rm=T))%>%
  arrange(-`Anteil komplett glücklich`)

#base R

ess$komplett <- ifelse(ess$happy==10,1,0)
ess <- ess[!is.na(ess$komplett),]
agghapp <- aggregate(ess$komplett,list(ess$Land),FUN=mean)
agghapp[order(agghapp$x),]

#Um die Verwendung von Subsetberechnungen kennenzulernen
#__________________________________
#5. Welches Land hat die grösste Einbusse im Durchschnittsglück gehabt#######
#von 2008 bis 2010?
ess[,.(
  Dif=mean(happy[year==2010],na.rm=T)-
    mean(happy[year==2008],na.rm=T)),
  by=Land][
    order(Dif)][!is.na(Dif)]

##alternativ: m.happy-m.happy[-.N] (eine Zeile vorher)
#bei sample mit nur 2008 und 2010 (geordnet)

#dplyr
rank<-ess%>%
  group_by(Land)%>%
  summarise(Dif=mean(happy[year==2010],na.rm=T)-
              mean(happy[year==2008],na.rm=T))%>%
  arrange(Dif)%>%
  filter(!is.na(Dif))

#Base R
zwei8 <- ess[ess$year==2008&!is.na(ess$happy),]
zwei10 <- ess[ess$year==2010&!is.na(ess$happy),]

agg8<-aggregate(zwei8$happy,list(zwei8$Land),FUN=mean)
names(agg8)<-c("Land","Happy.08")
agg10<-aggregate(zwei10$happy,list(zwei10$Land),FUN=mean)
names(agg10)<-c("Land","Happy.10")

ess.08.10 <- merge(agg8,agg10,by="Land",all.x = F) 
ess.08.10$dif <- ess.08.10$Happy.10-ess.08.10$Happy.08
ess.08.10[order(ess.08.10$dif),c("Land","dif")]


library(ggplot2)
ggplot(rank,aes(reorder(Land, Dif),Dif))+#einfache Art Faktoren nach Grösse einer anderen Dimension zu sortieren
  geom_bar(stat = "identity")+
  theme(axis.text.x = element_text(angle = 90, hjust = 1))+
  labs(x="",y="Glücklichkeitsveränderung 2010-2008")


# 6. Berechnen sie gleichzeitig den Modus bezüglich Glücklichkeit und Support für Umverteilung#####
# Um Funktion kennenzulernen, mit der man mehrere Variablen gleichzeitig modifiziert:
Mode <- function(x) {
  names(which.max(table(x)))
}

ess[,lapply(.SD, Mode),Land, .SDcols=c("happy","gincdif")]

#dplyr
ess %>%
  group_by(Land)%>%
  summarise_at(.vars = vars(happy, gincdif),
            .funs = Mode)

#Base R
aggregate(cbind(happy, gincdif)~Land,FUN = Mode,data = ess)

17.3 Vergleich Base-R, dplyr und data.table

# Beispiel data.frame:
data <- data.frame(A=1:10, B=2:11, C=3:12, D=c("a","b"))

# Beispiel data.table:
library(data.table)
data <- data.table(A=1:10, B=2:11, C=3:12, D=c("a","b"))
Base-R dplyr data.table
Auswahl von Zeilen nach Kriterium B==2 data <- data[data$B==2, ] data <- filter(data, B==2) oder
data <- data %>% filter(B==2)		 data <- data[B==2] oder
data <- data[B==2, ]
Auswahl von Spalten B und C data <- data[c("B", "C")] oder
data <- data[, 2:3] oder
data <- data[, c("B", "C")] oder
vars <- c("B", "C") data <- data[, vars]		 data <- select(data, B, C) oder
data <- data %>% select(B,C)		 data <- data[, .(B, C)] oder
data <- data[, c("B", "C")] oder
vars <- c("B", "C") data <- data[, ..vars]
Erstellen neuer Variablen E (A/B) data$E <- data$A/data$B data$E <- mutate(data, A/B) oder
data <- data %>% mutate(E = A/B)		 data[,E := A/B] oder
data[, E:=A/B, ] oder
data[, ":="(E=A/B)]
Erstellen neuer, gruppierter Variablen E (Mittelwert A gruppiert nach D) data$E <- ave(data$A,data$D,FUN=mean) data <- group_by(data,D)
data <- mutate(data,E = mean(A)) oder
data <- data %>%
   group_by(D) %>%
  mutate(E = mean(A))		 data[,E := mean(A),D] oder
data[,E := mean(A),by=D] oder
data[, ":="(E=mean(A)),by=D]
Aggregieren Mittelwert von A und B sapply(data[c("A", "B")], mean) data %>%
  summarise_at(vars=vars(A, B), .funs = mean)		 data[, .(meanB = mean(B), meanA = mean(A))] oder
data[, lapply(.SD, mean), , .SDcols = c("A", "B")]
Nach B gruppiert Aggregieren durch Mittelwert von A und C aggregate(data$A,data$B,mean) aggregate(data$C,data$B, mean) data <- group_by(data,B)
data <- summarise_at(data,vars = vars(A,C), .funs = mean) oder
data %>%
  group_by(B) %>%
  summarise_at(vars = vars(A,C), .funs = mean)		 data[, .(meanC = mean(C), meanD=mean(A)), B] oder
data[, lapply(.SD, mean),B, .SDcols=c("C", "A")]

17.4 Tidy vs Messy Data ↑

  • Hadley Wickham hat versucht, einige Konventionen zur sauberen Darstellung/Ablage von Daten aufzuschreiben
  • siehe http://vita.had.co.nz/papers/tidy-data.pdf und http://tidyverse.org/
  • Quintesenzen:
    • sobald Daten einmal sauber (tidy) sind, können Analysetools (Plotting, Modellfitting) auch sauber und ohne Zusatzaufwand arbeiten (z.B. ggplot2, reshape, lm/glm)
    • Fälle in Zeilen, Variablen (Eigenschaften dieser Fälle) in Spalten
  • Daten sind Messy, wenn
    • Spalten nicht/inkonsistent beschriftet sind
    • Eine Spalte mehr als eine Variable beinhaltet
    • Variablen statt in Spalten auch in Zeilen auftauchen
    • unterschiedliche Beobachtungseinheiten in derselben Tabelle sind
    • aber auch wenn wir den falschen Datentyp (z.B. String statt Datum, Character statt Numeric, etc.) verwenden.
  • wir brauchen also Werkzeuge um
      1. Strings und Werte zu manipulieren und
      1. Daten richtig anzuordnen

17.4.1 Tidy Data

  • Tidy Data liegen vor, wenn Beobachtungen in Zeilen sind und Attribute dieser Beobachtungen in Spalten
  • Sind die Daten anders organisiert, können wir melt() und dcast() aus dplyr resp. data.table verwenden
  • Es gibt neuere “Vokabeln” wie gather und spread aus tidyr sowie pivot_longer und pivot_wider ebenfalls aus tidyr
  • Alles macht im wesentlichen dasselbe
  • Installation des tidyr Packages macht generell Sinn. Dieses enthält eine Vielzahl von weiteren Paketen, die wir später teilweise noch brauchen werden (z.B. dplyr)

17.4.2 Beispiel: Tidy vs Messy Data ↑

  • Betrachten wir den Beispieldatensatz billboard. Dieser enthält Chartplatzierungen von Songs für 76 Wochen (jeweils in separaten Spalten).
  • Angenommen wir wollen die Wochen nicht nebeneinander (wide) sondern untereinander (long) haben?
  • Das aktuell beste Tool dafür ist pivot_longer (Gegenstück ist pivot_wider) aus Package tidyr.
library(tidyr)
billboard

# reshapen von wide nach long:
billboard %>% 
  pivot_longer(
    cols = starts_with("wk"), 
    names_to = "week", 
    values_to = "rank",
    values_drop_na = TRUE
  )

# diese Darstellung macht mehr Sinn, wenn wir Auswertungen der Chartplatzierungen machen wollen.

# Ein bisschen Finetuning:
billboard %>% 
  pivot_longer(
    cols = starts_with("wk"), 
    names_to = "week", 
    names_prefix = "wk",
    names_transform = as.integer,
    values_to = "rank",
    values_drop_na = TRUE,
  )
us_rent_income

# Daten sind nicht tidy weil estimate und moe Informationen zu verschiedenen Attributen (income und rent) mischen!

us_rent_income %>% 
  pivot_wider(
    names_from = variable, 
    values_from = c(estimate, moe)
  )

# beide Variablen werden aufgesplittet in Attribute bzgl income und rent.

17.4.3 Übung: Tidy data (interaktiv) ↑

  1. Laden Sie die Daten http://www.farys.org/daten/satisfaction.csv. Es handelt sich um erfundene Daten die das Alter (age) und die Zufriedenheit (c44) von 5 Personen über die Jahre 2012 bis 2014 enthalten. Sie sehen: es gibt für jedes Beobachtungsjahr und jede Beobachtungsdimension (Alter und Zufriedenheit) eine Variable. Die Information “Beobachtungsjahr” steckt implizit in den Variablennamen obwohl das Jahr eine eigene Variable sein sollte.
  2. Erstellen sie einen “tidy” Datensatz im long-Format. Dieser soll vier Spalten haben - ID, Jahr, Alter, Zufriedenheit - und eine Zeile für jede Messung einer Person in einem Jahr.

Das fertige Ergebnis soll so aussehen:

> data.wide
# A tibble: 15 x 4
   idpers  jahr zufriedenheit alter
    <int> <dbl>         <int> <int>
 1      1  2012             5    45
 2      1  2013             5    46
 3      1  2014             6    47
 4      2  2012             6    57
 5      2  2013             9    58
 6      2  2014             9    59

Tipps:

  • Erstellen sie zuerst einen “very long” Datensatz mit drei Variablen (idpers, name und value; z.b. mit pivot_longer).
  • Erstellen sie eine Variable “Dimension”, die “Alter” oder “Zufriedenheit” als Werte haben kann.
  • Hierfür können Sie z.B. ifelse plus grepl() verwenden oder auch andere Tools aus der Regex-Werkzeugkiste.
  • Erstellen Sie eine Variable “jahr”, z.B. durch Extrahieren der Jahresinformation mit substr() oder anderen Werkzeugen.
  • Überführen Sie zuletzt den Datensatz in ein tidy-Format (für jede Dimension eine Spalte, mit pivot_wider).

17.4.4 Lösung: Tidy data ↑

data <- read.csv("http://www.farys.org/daten/satisfaction.csv")
library(tidyr)

# Daten von wide nach long
data.long <- pivot_longer(data, cols = p12c44:age14)

# Auf welche Messdimension bezieht sich der Messwert?
data.long$dimension <- ifelse( # weise der Spalte dimension einen von zwei Werten zu, je nachdem ob folgende Bedingung erfüllt ist:
  grepl("c44", data.long$name), # Wenn die Spalte mit den ursprünglichen Spaltennamen ein «c44» enthält, 
  "zufriedenheit", # dann handelt es sich um eine Messung der Lebenszufriedenheit, also weise "zufriedenheit" zu, 
  "alter") # wenn nicht, dann weise "alter" zu.

# Auf welches Beobachtungsjahr bezieht sich der Messwert?
data.long$jahr <- ifelse(data.long$dimension == "zufriedenheit", # Hat die Spalte den Wert Lebenszufriedenheit
                    2000 + as.numeric(substr(data.long$name, 2, 3)), # dann ist die Jahresinformation and 2. Bis 3. Stelle
                    2000 + as.numeric(substr(data.long$name, 4, 5))) # sonst (Alter) an 4. Bis 5. Stelle

# Reshape ins wide Format damit es tidy wird:
data.long$name <- NULL # diese Spalte sollte entfernt werden, da für pivot_wide nicht klar ist, was es damit tun soll!
data.wide <- pivot_wider(data.long, names_from = dimension)



# Alternativ (und viel besser): falls die Namen in einem konsistenteren Format wären könnte man sich etwas Arbeit sparen:

# Bauen wir die Namen ein bisschen um:
names(data) <- sub(
  "p([0-9]{2})c44", 
  "satisfaction\\2\\1", 
  names(data))

# danach können wir die Struktur der Daten gut durch ein pattern beschreiben.

pivot_longer(data, cols = satisfaction12:age14, 
             names_pattern = "(.*)([0-9]{2})",
             names_to = c(".value", "year"))

# .value erklärt pivot_longer, dass dieser Teil den gemessenen Wert beschreibt und jeweils die neu entstehende Spaltenbeschriftung sein soll.

17.4.5 Beispiel: Von Messy zu sauber mit diversen Tools ↑

# Wetterdaten
weather <- read.table("https://raw.githubusercontent.com/justmarkham/tidy-data/master/data/weather.txt", header=TRUE)
head(weather) # hier sind Variablen in Zeilen und Spalten

# Daten reshapen von wide nach long
library(tidyr)
weather1 <- pivot_longer(weather, cols = d1:d31, values_drop_na = TRUE)
head(weather1)

# saubere Spalte für "day"
library(stringr)    # für str_replace(), str_sub()
weather1$day <- as.integer(str_replace(weather1$name, "d", ""))

# die krude Spalte "variable" brauchen wir nicht. Löschen durch Zuweisen von NULL.
weather1$name <- NULL

# die Spalte element beherbergt zwei unterschiedliche Variablen tmin und tmax. 
# Diese sollen in zwei Spalten:
weather1$element <- tolower(weather1$element) # Kleinbuchstaben
weather.tidy <- pivot_wider(weather1, names_from = element)
head(weather.tidy)

# das Datum laesst sich zudem in einer Spalte darstellen als echtes Datum:
weather.tidy$date <- as.Date(paste(weather.tidy$year, 
                                   weather.tidy$month, 
                                   weather.tidy$day, sep="-"))
weather.tidy[, c("year", "month", "day")] <- NULL
head(weather.tidy)


# ABER: was für ein Mess mit weather1 weather.tidy etc. Warum nicht Piping verwenden?

weather %>%
  pivot_longer(cols = d1:d31, values_drop_na = TRUE) %>%
  mutate(day = as.integer(str_replace(name, "d", "")),
         element = tolower(element)) %>%
  select(-name) %>%
  pivot_wider(names_from = element) %>%
  mutate(date = as.Date(paste(year, month, day, sep = "-"))) %>%
  select(-year, -month, -day)

17.5 Mergen von Daten ↑

  • Wird in der Praxis sehr häufig gebraucht. Klassisch gibt es merge() (siehe Quick R)
  • Die Join Funktionen in dplyr sind allerdings schneller
  • inner_join(x, y, by = NULL, copy = FALSE, ...) # behält alle Fälle die zwischen den Daten x und y überlappen; häufigster Fall
  • left_join(x, y, by = NULL, copy = FALSE, ...) # behält alle Fälle von x. Missing falls in y nicht vorhanden
  • semi_join(x, y, by = NULL, copy = FALSE, ...) # behält alle Fälle von x und behält nur die Spalten von x
  • anti_join(x, y, by = NULL, copy = FALSE, ...) # alle Fälle, für die es keinen Match gibt
    • x = Datensatz 1, y = Datensatz 2
    • by = “Matchingvariable (z.B. ID, Name, Countrycodes, etc.)”
    • ?join: Das Hilfe File ist sehr zu empfehlen!
  • ebenfalls häufig gebraucht: rbind() um data frames mit gleichen Spalten zeilenweise zusammenzukleben

17.5.1 Beispiel: Mergen von Daten ↑

name <- c("Rudi","Simon","Daniela","Viktor")
geschlecht <- c("Mann", "Mann", "Frau","Mann")
daten1 <- data.frame(name, geschlecht)

name <- c("Johanna","Rudi","Simon","Daniela")
alter <- c(33,32,38,45)
daten2 <- data.frame(name, alter)

daten1
daten2

# klassisch
merge(daten1, daten2, by="name") # entspricht inner_join)

# oder über match()
# daten1$alter <- daten2$alter[match(daten1$name,daten2$name)] # entspricht left_join

# mit dplyr
inner_join(daten1,daten2,by="name")
left_join(daten1,daten2,by="name")
right_join(daten1,daten2,by="name")
full_join(daten1,daten2,by="name")

# mit data.table
# mit data.table
name <- c("Rudi","Simon","Daniela","Viktor")
geschlecht <- c("Mann", "Mann", "Frau","Mann")
daten1 <- data.table(name, geschlecht)

name <- c("Johanna","Rudi","Simon","Daniela")
alter <- c(33,32,38,45)
daten2 <- data.table(name, alter)
setkey(daten1, name)
setkey(daten2, name)

#Inner join
daten1[daten2,nomatch=0]

#Right-join
daten1[daten2]

#Left-join
daten2[daten1]

#Full outer-join
unique_keys <- unique(c(daten1[, name], daten2[, name]))
daten1[daten2[.(unique_keys)]]
#data.table kennt auch eine verbesserte Version von merge
merge(daten1, daten2, all = T)

# rbind / zeilenweise verknüpfen von Daten
name <- c("Rudi","Simon","Daniela","Viktor")
geschlecht <- c("Mann", "Mann", "Frau","Mann")
daten1 <- data.frame(name, geschlecht)

name <- c("Johanna","Ralf")
geschlecht <- c("Frau","Mann")
daten2 <- data.frame(name, geschlecht)

rbind(daten1,daten2)

# rbind() ist sehr langsam. Besser ist rbindlist() aus data.table
rbindlist(list(daten1,daten2))

18 Datenimport und -export ↑

18.1 Häufige Packages und Funktionen zum Datenimport und -export ↑

  • read.table(): Liest Textformate (.txt, .csv)
    • read.table("meinedaten.csv", header = TRUE, sep = ",", dec = ".")
      • header = TRUE: wenn die erste Zeile Spaltennamen enthält
      • sep = ",": Textseparator
      • dec = ".": Dezimalzeichen (Default ist .)
      • Für mehr siehe die Hilfe ?read.table
  • Text (.csv oder .txt) ist der kleinste Gemeinsame Nenner zwischen Softwarepaketen und empfiehlt sich, wenn mit vielen verschiedenen Programmen gearbeitet wird
  • Für sehr grosse Dateien empfiehlt sich die Funktion fread() aus data.table (Geschwindigkeit ca. Faktor 10)
  • Oft liegen andere Formate vor weil primär eine bestimmte Software eingesetzt wird, z.B. STATA, SPSS, SAS
    • foreign: kann enbtsprechende Formate lesen, z.B. (.dta, .por)
    • haven: neueres Package, das versucht foreign abzulösen
  • Ebenfalls verbreitet ist xls/xlsx. Hier gibt es mehrere Möglichkeiten, die sinnigsten sind:
    • read_excel() aus Package readxl: baut auf C. Ebenfalls sehr unproblematisch und versteht exceltypische Sheet-/Zellennotation.
    • read.xlsx() aus Package openxlsx (keine Abhängigkeiten jenseits von R, liest aber keine alten .xls Files, nur xlsx)
    • read.xls() aus Package gdata, braucht allerdings eine Perl Installation
    • read.xlsx() und read.xlsx2() aus Package xlsx; kann auch bestimmte Zeilen-/Spaltenbereiche auslesen statt komplette Sheets, benötigt jedoch Java; ggf. Java installieren und sicherstellen, dass R und Java dieselbe Architektur (32 oder 64bit) haben. Für Probleme bei Mac-Usern: http://stackoverflow.com/questions/28796108/r-rstudio-yosemite-and-java
  • Daten exportieren
    • Idealerweise als Text
    • write.table(swiss, "swiss.txt", sep=",") speichert das Objekt swiss kommasepariert nach swiss.txt
    • Siehe ?write.table für weitere Optionen (z.B. Zeilennamen, Separator, Dezimalzeichen)
    • Für eine gute Übersicht über andere Formate siehe Quick R

18.2 Übung: Datenimport und -export / Mergen / Tidy Data ↑

  1. Beim BfS gibt es Daten zur ständigen Wohnbevölkerung der Schweizer Kantone. Lesen Sie Daten über die Population (Spalte B) der 26 Schweizer Kantone (Spalte A) ein. Diese Befinden sich in Sheet “2022” des Excel-Sheets “Bilanz der ständigen Wohnbevölkerung nach Kantonen” (ggf. nur relevante Zeilen auswählen!). Hinweis: Einlesen direkt via URL ist nicht möglich, da xlsx kein Plaintext ist sondern ein “binary zip”. Laden Sie die Datei entweder manuell herunter oder (besser) mit Hilfe von download.file(..., mode="wb") (write binary mode). Hinweis 2: Achten Sie auf das File-Encoding der Daten.
  2. Gehen Sie wie in 1. vor und besorgen Sie von der Website des BfS zum Thema Forstwirtschaft das Excel Sheet “Waldflächen nach Kantonen” und lesen Sie Kantonsnamen und die Gesamtwaldfläche in Hektar ein. Hinweis: Bei den Daten wurde eine konfuse Anzahl Leerzeichen hinter die Kantonsnamen gehängt (read_excel trimmt diese automatisch weg, read.xlsx jedoch nicht). Säubern Sie die Daten, damit die Kantonsnamen mit denen aus 1. korrespondiert.
  3. Verknüpfen (mergen) Sie die beiden Datensätze über die Kantonsnamen. Können alle Kantone gematcht werden?
  4. Berechnen Sie für alle Kantone die Anzahl Bäume pro Person (1 Hektar = 400 Bäume) und geben Sie eine Rangliste aus (z.B. mit arrange() aus dplyr).
  5. Speichern Sie diese Rangliste als csv (z.B. mit write.table()).

Tipps: Schauen Sie in die Hilfe von read_excel(). Hilfreiche Argumente, die hier spezifiziert werden können sind z.B. range, sheet und col_names.

18.3 Lösung: Datenimport und -export / Mergen / Tidy Data ↑

# 1.
library(readxl)

# Datei runterladen als binary
download.file("https://dam-api.bfs.admin.ch/hub/api/dam/assets/26565330/master", 
              "bevoelkerung.xlsx", 
              mode="wb")

# einlesen als UTF-8 kodiertes File. Sheet 2014, Zeile 8-33, 
# ohne Header (den setzen wir besser manuell). über colIndex 
# nehmen wir nur die Spalten 1 (A) und 2 (B)


bevoelkerung <- read_excel("bevoelkerung.xlsx", sheet = "2022", range = "A6:B31", col_names = c("Kanton", "Einwohner"))


# 2.
download.file("https://dam-api.bfs.admin.ch/hub/api/dam/assets/26305060/master",
              "wald.xlsx", 
              mode="wb")

wald <- read_excel("wald.xlsx", sheet = 1, range = "A6:C35", col_names = c("Kanton", "Eigentuerm", "Waldflaeche"))[, c(1, 3)]


# säubern (entfällt bei readxl im Vergleich zu read.xlsx)
library(stringr)
wald$Kanton <- str_trim(wald$Kanton, side="both")

# 3. / 4.
# z.B. über dplyr
library(dplyr)
left_join(bevoelkerung, wald, by = "Kanton")

# St.Gallen und die beiden Appenzells werden nicht gematcht 
# da die Schreibweise minimal anders ist, fix:

wald$Kanton <- str_replace(wald$Kanton, pattern = "\\. ", replacement = ".") 

# der Punkt wird hier escaped, da es sich normal um einen
# Platzhalter für reguläre Ausdrücke handelt
# besser ist normalerweise über das Kantonskürzel oder
# Nummer zu matchen, dann stört die unterschiedliche Schreibweise nicht.

waldranking <- left_join(bevoelkerung,wald,by="Kanton") %>%
  mutate(bpp = Waldflaeche*400/Einwohner) %>%
  arrange(-bpp)

waldranking

# 5. 
write.csv(waldranking, "waldranking.csv")

18.4 Daten aus Datenbanken lesen

  • relevante Packages: RODBC, RJDBC, RMySQL, RPostgreSQL, DBI
  • ein kleines Beispiel:
library(DBI)
library(RMySQL)

# Verbindung definieren
con <-  dbConnect(MySQL(), 
                  username = "dataviz", 
                  password = "CASdataviz2016",
                  host = "db4free.net", 
                  port = 3306,
                  dbname = "cas_dataviz"
)

# Get-Anfrage schicken
dbGetQuery(con, "show databases")

# neue Tabelle erstellen
dbSendQuery(con, "CREATE TABLE anmeldungen (id INT(6) UNSIGNED AUTO_INCREMENT PRIMARY KEY, 
  vorname VARCHAR(30) NOT NULL,
  nachname VARCHAR(30) NOT NULL,
  email VARCHAR(50),
  reg_date TIMESTAMP)"
) 

# schauen welche Tabellen es gibt
dbGetQuery(con, "show tables")

#   Tables_in_cas_dataviz
# 1           anmeldungen

# Alle Zeilen der Tabelle löschen
# dbSendQuery(con, "DELETE from anmeldungen")

dbSendQuery(con, "INSERT into anmeldungen (vorname, nachname) values ('Rudi', 'Farys')")
dbSendQuery(con, "INSERT into anmeldungen (vorname, nachname,email) values ('Test', 'Test','test@test.de')")

# einfache Abfrage
daten <- data.frame(dbGetQuery(con, "SELECT * FROM anmeldungen"))

daten

dbDisconnect(con)

19 Webscraping ↑

19.1 HTML/XML data ↑

19.1.1 Beispiel: HTML Tabelle einlesen ↑

library(rvest)
raw <- read_html("https://en.wikipedia.org/wiki/Cantons_of_Switzerland")
table <- html_table(raw, fill = TRUE, header = TRUE)[[3]]
head(table)

19.2 CSS Selectors für komplexere HTML-Strukturen ↑

  • CSS Selektoren sind eine Syntax zur Definition und Navigation durch HTML/XML-Dokumente.
  • Bezugnahme auf Elemente (z.B. einer Website oder eines XML-Dokuments) durch deren Pfad, Position oder Bedingungen.
  • Wenn Sie Informationen von einer Website extrahieren: Element markieren, mit der rechten Maustaste klicken und “Inspizieren” / “Element prüfen” auswählen oder ein entsprechendes Browser-Plugin verwenden.
  • Paket: am besten rvest, das eine Reihe von Funktionen zum komfortablen Extrahieren von Elementen enthält, z.B. html_table, html_node``,html_attr`` und andere.
  • Das Lesen der Website in R erfolgt in der Regel durch read_html()``` vonxml2, auf dem ``rvest aufbaut.

19.2.1 Beispiel: CSS Selectors ↑

library(rvest)
# Take a look at
# https://www.transfermarkt.com/fill/profil/spieler/28003 and
# https://www.transfermarkt.com/fill/profil/spieler/401173

# Notice how its possible to navigate through player profiles when we use a list of player IDs, e.g.
listofplayers <- c(28003, 401173)

for(playerid in listofplayers) {
  url <- paste0("https://www.transfermarkt.com/fill/profil/spieler/", playerid)
  raw <- read_html(url)
  ...
}

# or even more performant with a function
# lets say we want to extract name and market value of the player.

# inspecting the element of market value reveals, that we have an anchor (<a>) within a <div> container that has class "data-header__market-value-wrapper" which is probably usable. The text within <a> is what we want:


raw <- read_html("https://www.transfermarkt.de/fill/profil/spieler/28003")
rawelement <- html_node(raw, ".data-header__market-value-wrapper") %>%
  html_text()
rawelement

# the dot means "class", but we can also be more explicit:
html_node(raw, "[class = 'data-header__market-value-wrapper']") %>% html_text()

# the raw text we need to tidy up, e.g.:
rawelement %>% 
  strsplit(., ",", fixed = TRUE) %>%
  .[[1]] %>% .[1] %>%
  gsub("€","",.) %>%
  as.numeric

# Players name:
# there are many many ways but e.g. open source code and look for players name. Select any suitable entry that you can easily tackle with css selectors. E.g. the title tag:

raw %>% 
  html_node("title") %>%
  html_text %>%
  strsplit(., " - ") %>%
  .[[1]] %>% .[1]

# in addition to extracting text using html_text we can also extract attributes, which is often helpful. Consider the following part of the source code:
# <meta property="og:title" content="Lionel Messi - Player profile 21/22" />

# we can use html_attr() here to have an alternative for extracting the players name from the Open Graph tags in the head of the html (see http://ogp.me/)
raw %>% html_node("[property='og:title']") %>% html_attr("content")

# we can now pack everything we want into a function and apply to lots of websites!

getPlayerStats <- function(playerid) {
  
  url <- paste0("https://www.transfermarkt.com/fill/profil/spieler/", playerid)
  raw <- read_html(url)
  
  marketvalue <- raw %>% 
    html_node(".data-header__market-value-wrapper") %>%
    html_text() %>% 
    strsplit(., ".", fixed = TRUE) %>%
    .[[1]] %>% .[1] %>%
    gsub("€","",.) %>%
    as.numeric
  
  name <- raw %>% 
    html_node("title") %>%
    html_text %>%
    strsplit(., " - ") %>%
    .[[1]] %>% .[1]
  
  data.frame("id" = playerid, "name" = name, "marketvalue" = marketvalue)
}

# Messi example
getPlayerStats(28003)

# now for a list of IDs
do.call(rbind, lapply(listofplayers, getPlayerStats))

19.2.2 Übung: CSS Selectors ↑

  • Werfen Sie einen Blick auf die Website https://www.ariva.de/DE0007100000
  • Wie Sie sehen können, können Sie als Teil der URL die ISIN (International Securities Identification Number) eingeben, um die Aktienkurse und andere Informationen über eine Aktie zu erhalten.
  • Schreiben Sie eine Funktion, die den aktuellen Kurs einer Aktie über ihre ISIN extrahiert.
  • Hinweise: Schauen Sie sich den HTML-Quellcode an, suchen Sie nach dem Preis und analysieren Sie die Struktur des Elements. Welche CSS-Selektoren können Sie anwenden, um den Preis zu extrahieren?
  • Bonus: Stellen Sie auch sicher, dass die Ausgabe Ihrer Funktion numerisch ist.

19.2.3 Lösung: CSS Selectors ↑

library(dplyr)
library(rvest)

arivaQuote <- function(isin) {
  url <- paste0("https://ariva.de/", isin)
  qt <- read_html(url) %>%
    html_node("[itemprop='price']") %>% 
    html_attr("content") %>% # alternativ ggf. auch html_text() plus ein bisschen trimming/säubern
    as.numeric()
  qt / 100
}

arivaQuote("DE0005190003") # BMW

# auch noch lustig: obwohl die Funktion nur 1 URL aufs mal verarbeiten kann, könnten wir trotzdem eine vektorisierte Version bauen:

arivaQuote <- Vectorize(arivaQuote)
isins <- c("DE0005190003", "US88160R1014") # BMW und Tesla

mydf <- data.frame(isins)
mydf$quote <- arivaQuote(mydf$isins)

19.3 Simulation eines Webbenutzers für sehr haarige Fälle↑

  • Webscraping wird von Tag zu Tag schwieriger, hier sind einige Gründe:
    • viele Websites verwenden AJAX, welches Daten erst dann in das Layout der Website lädt, nachdem der Quellcode geladen wurde. Das Lesen des Quellcodes ergibt daher nur eine leere Struktur einer Website ohne Daten.
    • Möglicherweise müssen Sie sich bei einer Website anmelden.
    • Datenbesitzer wollen ihre Daten ggf. monetarisieren und wollen daher nicht, dass gescrapt wird.
  • RSelenium ist eine Lösung für einige dieser Probleme, da man damit einen Browser-Benutzer simulieren kann. Auf diese Weise können Sie z.B:
    • eine Login-Seite überwinden und Cookies, Session-IDs, die benötigt werden, verwenden.
    • Daten laden, die nicht im reinen Quellcode vorhanden sind (z.B. wenn Daten beim Scrollen dynamisch geladen werden).
  • Bevor Sie Websites scrapen, die schwer zu scrapen sind, kontaktieren Sie am besten den Website-Besitzer und bitten diesen um Unterstützung. Oftmals bekommt man dieselben Daten dadurch bereits in maschinenlesbarer Form.
  • RSelenium zu zeigen, geht über den Rahmen dieses Kurses hinaus, aber wenn nötig, können wir Ihnen Beispielcode zur Verfügung stellen.

20 APIs - Application programming interface ↑

20.1 Beispiel: World Bank Indicators

20.2 JSON-Daten ↑

  • JSON (java script object notation) ist das populärste Datenformat für den Datenaustausch über APIs
  • In R gibt es mehrere Packages hierfür. Das einschlägigste Package ist jsonlite.
library(jsonlite)

# Mini Beispiel aus dem Package
json <-
'[
  {"Name" : "Mario", "Age" : 32, "Occupation" : "Plumber"}, 
  {"Name" : "Peach", "Age" : 21, "Occupation" : "Princess"},
  {},
  {"Name" : "Bowser", "Occupation" : "Koopa"}
]'
mydf <- fromJSON(json)
mydf


# editieren wir die Daten ein bisschen
mydf$Ranking <- c(3, 1, 2, 4)

# und wandeln es zurück nach JSON
toJSON(mydf, pretty=TRUE)
  • Wie könnte das ganze für das Weltbank Beispiel aussehen?
# Indicator: NY.GDP.PCAP.PP.CD (GDP per capita PPP)
url <- "http://api.worldbank.org/v2/country/all/indicator/NY.GDP.PCAP.PP.CD?date=2021&format=json&per_page=10000"

# was für Teile gibt es hier in der URL?

baseurl <- "http://api.worldbank.org"
endpoint <- "/v2/country/all/indicator/NY.GDP.PCAP.PP.CD"
parameters <- "?date=2021&format=json&per_page=10000"
url <- paste0(baseurl, endpoint, parameters)

library(jsonlite)
res <- fromJSON(url)

str(res)
class(res[[2]])
names(res[[2]])

library(dplyr)
res[[2]] %>%
  select(country, gdp = value) %>%
  arrange(-gdp) %>%
  head(10)

20.3 XML-Daten ↑

  • Grundsätzlich das gleiche Konzept wie JSON, jedoch mit unterschiedlichem Markup.
# XML mini Beispiel
library(XML)
download.file("https://www.w3schools.com/xml/simple.xml", "simple.xml")
# Leider erfordern die meisten Seiten, dass Sie sie zuerst herunterladen, da xmlParse() https nicht unterstützt.
roh <- xmlParse("simple.xml")
liste <- xmlToList(roh)
df <- xmlToDataFrame(roh)

20.4 APIs, die bereits in einem R-Paket verpackt sind ↑

  • Einige Pakete wurden bereits für R geschrieben, die es einfacher machen, mit einzelnen APIs zu arbeiten.

  • Für die eben genannten World development indicators (Weltbank) gibt es ein Package WDI, aber auch für Daten der OECD, IWF-Daten, Börsendaten (quantmod und andere), twitteR, rfacebook, quandl, tidycensus, gapminder, Rspotify uvvm.

# install.packages("WDI") 
library(WDI)
WDIsearch(string="gdp", field="name", cache=NULL)
DF <- WDI(country="all", 
          indicator="NY.GDP.PCAP.PP.CD", # Pro-Kopf-BIP kaufkraftbereinigt
          start=2010, end=2015) 

# Welche sind die reichsten Länder im Jahr 2010?
library(dplyr)
filter(DF, year==2010) %>%
arrange(-NY.GDP.PCAP.PP.CD) %>%
  head()

# Und wie sieht es 2015 aus?
filter(DF, year==2015) %>%
arrange(-NY.GDP.PCAP.PP.CD) %>%
  head()

# Was passiert da im Hintergrund?

# ?WDI
# ?wdi.query

country = "all"
indicator = "NY.GDP.PCAP.PP.CD"
extra = FALSE
cache = NULL
latest = NULL
language = "en"
years <- NULL


url <- paste0("https://api.worldbank.org/v2/", language, 
       "/country/", country, "/indicator/", indicator, "?format=json", 
       years, "&per_page=32500", "&page=", 1:10, latest)

url

20.5 Übung APIs (optional für Interessierte)

  • Schauen Sie sich die Dokumentation zu SWAPI (The Star Wars API) an: https://swapi.dev/documentation
  • Machen Sie eine API-Abfrage (json-Format) um Daten für alle Figuren (Endpoint /people) herunter zu laden.
  • Nehmen Sie aus dem erhaltenen Datenobjekt den Namen und das Geburtsjahr und bauen Sie daraus einen data.frame.
  • Bonus: Sie erhalten nur die ersten 10 Treffer aber es gibt 82 Figuren. Wie könnte man einen data.frame für alle 82 erzeugen? (Hinweis: es gibt für die API einen Parameter page)
  • Bonus 2: Wie könnte man den Heimatplaneten (homeworld) mit in die Daten bekommen?

20.6 Lösung APIs

library(jsonlite)

starwars <- fromJSON("https://swapi.dev/api/people/?format=json")
df <- data.frame(name = starwars$results$name, geburtsjahr = starwars$results$birth_year)

# Für alle 82 brauchen wir 9 Abfragen:
df <- NULL
for(p in 1:9) {
  starwars <- fromJSON(paste0("https://swapi.dev/api/people/?format=json&page=", p))
  df <- rbind(df, data.frame(name = starwars$results$name, geburtsjahr = starwars$results$birth_year))
}

# Für den Heimatplanet brauchen wir jeweils eine eigene Abfrage:
df <- NULL
for(p in 1:9) {
  starwars <- fromJSON(paste0("https://swapi.dev/api/people/?format=json&page=", p))
  df <- rbind(df, data.frame(name = starwars$results$name, 
                             geburtsjahr = starwars$results$birth_year, 
                             heimatplanet = starwars$results$homeworld[1]))
}

# man könnte dafür z.B. eine eigene Funktion schreiben:
getHomeworld <- Vectorize(function(api_call) {
  res <- fromJSON(api_call)
  res$name
})

df$heimatplanet.name <- getHomeworld(df$heimatplanet)
# Achtung: das ginge auch effizienter (ohne API Calls doppelt zu machen)

21 Export von Tabellen ↑

21.1 Beispiel: Export von Tabellen

library(huxtable) # ggf. Installation von officer und flextable erforderlich
library(dplyr)
library(tidyr)

# Beispiel für Formatierung

mtcars %>%
  group_by(Zylinder = cyl, Gaenge = gear) %>%
  summarise(PS = round(mean(hp)),
            "Sek. 1/4 Meile" = round(mean(qsec))) %>%
  arrange(Zylinder, Gaenge) %>%
  as_hux(add_colnames = TRUE) %>%
  set_bold(1, 1:4, TRUE) %>%
  set_bottom_border(1, 1:4, TRUE) %>%
  set_left_border(everywhere,3, TRUE) %>%
  set_bold(1:9, 1:2 ,TRUE) %>%
  merge_cells(2:4, 1) %>%
  merge_cells(5:7, 1) %>%
  merge_cells(8:9, 1) %>%
  map_text_color(everywhere, 3:4, 
                 by_quantiles(c(0.2, 0.8), c("blue", "darkgreen", "red"))) %>%
  quick_docx(., file = "Test_format.docx")


# Beispiel für Summary Tabelle
tabelle <- read.csv("http://farys.org/daten/Prestige.csv") %>%
  select(education, income) %>%
  #Kann beliebig ergänzt werden
  summarise_all(list(Minimum = min,
                     Median = median,
                     Mittelwert = mean,
                     Maximum = max,
                     Fallzahl = ~sum(!is.na(.))), na.rm = TRUE) %>%
  pivot_longer(everything()) %>%
  mutate(Dimension = gsub("_.*", "", name),
         Mass = gsub(".*_", "", name),
         value = as.character(round(value))) %>%
  select(-name) %>%
  pivot_wider(names_from = Mass, values_from = value) %>%
  as_hux(add_colnames = T) %>%
  set_bold(1,everywhere,TRUE)

# Schöne Spaltennamen
tabelle[, 1] <- c("", "Bildung (in Jahren)", "Einkommen")
quick_docx(tabelle, file = "Test_summary.docx")



# Regressionstabellen

x1 <- rnorm(100)
x2 <- rnorm(100)
y <- rnorm(100) + 2*x1 + 1*x2

fit1 <- lm(y~x1)
fit2 <- lm(y~x1+x2)

huxreg(fit1, fit2)

#https://cran.r-project.org/web/packages/huxtable/vignettes/huxreg.html für Formatierungsmöglichkeiten

21.2 Übung: Export von Tabellen↑

  1. Öffnen Sie die Daten http://www.farys.org/daten/bmi.dta. Diese enthalten (unter anderem) Alter, Geschlecht und den BMI (Body-Mass-Index) aus einer Befragung.

  2. Erzeugen Sie eine Summary-Statistik für die Variablen BMI, Alter und Einkommen. Nutzen Sie am besten das Beispiel als Template und passen Sie dieses an.

  3. Hübschen Sie die Tabelle noch etwas auf, indem Sie die Zeilen besser beschriften.

  4. Exportieren Sie die Tabelle als Word/HTML, Latex oder Ascii

  5. Erzeugen/exportieren Sie eine Tabelle, die alle Variablen der fünf höchsten und der fünf niedrigsten BMI-Werte zeigt.

21.3 Lösung: Export von Tabellen ↑

# 1.
library(foreign)
bmi <- read.dta("http://www.farys.org/daten/bmi.dta")

# 2./3./4.
library(dplyr)

tabelle <- bmi %>%
  select(bmi, alter, einkommen)%>%
  #Kann beliebig ergänzt werden
  summarise_all(list(Minimum = min,
                     Median = median,
                     Mittelwert = mean,
                     Maximum = max,
                     Fallzahl = ~sum(!is.na(.))), na.rm = TRUE) %>%
  pivot_longer(everything()) %>%
  mutate(Dimension = gsub("_.*", "", name),
         Mass = gsub(".*_", "", name),
         value = as.character(round(value))) %>%
  select(-name) %>%
  pivot_wider(names_from = Mass, values_from = value) %>%
  as_hux(add_colnames = T) %>%
  set_bold(1,everywhere,TRUE)

#Schöne Namensgebung der Spalten
tabelle[, 1] <- c("Dimension", "Body-Mass-Index", "Alter", "Einkommen")
tabelle

quick_docx(tabelle, file = "summarytab.docx")


# 5. Es gibt viele Wege, Vorschlag:
topbottom <- rbind(arrange(bmi, -bmi) %>% head(n = 5), 
      arrange(bmi, bmi) %>% head(n = 5))

tabelle <- topbottom %>%
  as_hux(add_colnames = T ) %>%
  set_bold(1, 1:ncol(topbottom), TRUE) %>%
  set_bottom_border(1, 1:ncol(topbottom), TRUE)
tabelle[1, ] <- c("Body-Mass-Index","Alter","Geschlecht","Einkommen","Bildungsjahre","ID")
tabelle

# siehe ggf. auch https://hughjonesd.github.io/huxtable/reference/number_format.html für Number formatting

22 Exkurs: Workflow ↑

23 Lineare Regression ↑

23.1 Beispiel: Lineare Regression ↑

  • Daten über das Prestige von Berufen (Untersuchungseinheiten sind Berufe, nicht Individuen)
i prestige \(= y_{i}\) education \(= X_{i2}\) income \(= X_{i3}\)
1 gov.administrators 68.8 13.11 12351
2 general.managers 69.1 12.26 25879
3 accountants 63.4 12.77 9271
4 purchasing.officers 56.8 11.42 8865
5 chemists 73.5 14.62 8403
6 physicists 77.6 15.64 11030
7 biologists 72.6 15.09 8258
8 architects 78.1 15.44 14163
9 civil.engineers 73.1 14.52 11377
10 computer.programers 53.8 13.83 8425
11 economists 62.2 14.44 8049
… … … … …
  • Abhängige Variable: Pineo-Porter-Prestige für den jeweiligen Beruf (Werte stammen aus Befragungen)
  • Forschungsfrage: Haben Einkommen und Bildung einen Einfluss auf das Berufsprestige?
  • Erklärende Variablen
    • Education: durchschnittliche Anzahl Bildungsjahre der Leute im jeweiligen Berufsfeld
    • Income: durchschnittlicher Lohn im Berufsfeld
  • Datensatz heisst Prestige und befindet sich im Package car (siehe ?Prestige) (“Companion to Applied Regression” - John Fox)
# Berufliches Prestige und Bildung bzw. Einkommen

# Paket "car" laden bzw. installieren, da dort der 
# Beispieldatensatz "Prestige" enthalten ist

install.packages("car")
library(car) 

# alternativ, falls car Probleme mit Abhaenigkeiten macht: 
# Prestige <- read.csv("http://farys.org/daten/Prestige.csv")

# Für folgendes Beispiel brauchen wir jedoch die Funktion scatter3d() aus car. 
# Das Dependency Problem lässt sich ggf. so lösen:
# install.packages("lme4") # dependency für altes pbkrtest
# packageurl <- "https://cran.r-project.org/src/contrib/Archive/pbkrtest/pbkrtest_0.4-4.tar.gz"
# install.packages(packageurl, repos=NULL, type="source") # von hand installieren
# install.packages("car") # jetzt car installieren mit den 
# manuell installierten dependencies
# library(car)

# Wie kann man sich eine Regression mit zwei erklärenden Variablen vorstellen?
# Als Ebene durch eine 3d Punktewolke!

scatter3d(Prestige$income,Prestige$prestige,Prestige$education, fit="linear") 
  
# Ein kleines Modell schätzen:
fit <- lm(prestige ~ education + income, data=Prestige)

# summary() ist eine generische Funktion. Für ein lm() Objekt wird 
# ein typischer Überblick über Koeffizienten und Modellgüte gegeben
  
summary(fit)

# Was steckt im fit Objekt?
names(fit)

# Zugriff auf Koeffizienten und andere Bestandteile
fit$fitted.values # Vorhersage (y-Dach)
fit$residuals     # Residuen = Fehler, d.h. die Abweichung der beobachteten Werte 
                    # von der Modellvorhersage
coef(fit)         # Koeffizienten
fit$coefficients

24 Andere wichtige Packages/Funktionen für statistische Modelle und mehr

24.1 OLS Diagnostik - Annahmen testen ↑

  • Ausreisseranalyse: avPlots() aus Package car
  • Linearität: component plus residual plot: crPlots() aus car
  • Unabhängigkeit der Residuen (Autokorrelation):
    • (partielle) Autokorrelationsfunktion: acf(), pacf()
    • Q-Test/Ljung-Box-Test: Box.test()
  • Homoskedastizität
    • ncv.test(fit) und spread.level.plot(fit) aus car
    • robuste Standardfehler lassen sich mit Funktionen aus dem Package sandwich berechnen.

24.2 Zeitreihenanalyse

# install.packages("quantmod")
library(quantmod)
getSymbols("^SSMI")
head(SSMI)

# install.packages("forecast")
library(forecast)

fit <- auto.arima(to.monthly(SSMI)[,"SSMI.Adjusted"])

plot(forecast(fit,h=6))

24.3 Data Mining

25 Take-Home-Exam: ↑

25.1 Ebay-Auktionen

  1. Laden Sie die Daten http://www.farys.org/daten/ebay.dta. Es handelt sich um Ebaydaten von Mobiltelefonauktionen. Die Variablen sind:
  • sold: Ob das Mobiltelefon verkauft wurde
  • price: Der erzielte Verkauftspreis
  • sprice: Der Startpreis der Auktion
  • sepos: Anzahl positiver Bewertungen des Verkäufers
  • seneg: Anzahl negativer Bewertungen des Verkäufers
  • subcat: Das Modell des Mobiltelefons
  • listpic: Kategorialer Indikator, ob die Auktion ein Thumbnail, ein “has-picture-icon” oder kein Thumbnail hat.
  • listbold: Dummy, ob die Auktion fettgedruckt gelistet ist
  • sehasme: Dummy, ob der Verkäufer eine “Me-page” hat oder nicht

  1. Sie interessieren sich dafür, ob Mobiltelefone bei Auktionen einen höheren Preis erzielen, wenn der Verkäufer des Geräts eine gute Bewertung hat. Erzeugen Sie eine Variable \(rating\), die den Anteil positiver Bewertungen an den Gesamtbewertungen misst. Schliessen Sie Fälle aus den Daten aus, für die weniger als 12 positive Bewertungen vorliegen.

  2. Bilden Sie zudem eine Variable \(makellos\) (TRUE/FALSE), ob der Verkäufer ein makelloses Rating (\(>98\%\) positive Bewertungen) hat oder nicht.

  3. Zeichnen Sie einen farblich geschichteten Boxplot: Y-Achse=Preis, X-Achse=Gerätetyp, farblich geschichtet nach Bewertung (makellos=grün sonst=rot). Orientieren Sie sich z.B. am Vorschlag von “Roger Bivand” im Helpfile zu boxplot(). Exportieren Sie die Grafik als PDF. Erzielen (rein optisch) die Verkäufer mit makellosem Rating einen höheren Verkaufspreis?

  4. Rechnen Sie zwei kleine Regressionsmodelle für den Preis von verkauften Geräten. Modell 1 soll als Prädiktoren den Modelltyp und das Rating beinhalten. Modell 2 soll zusätzlich die Variable listpic beinhalten. Haben das Rating und die Thumbnails einen Einfluss auf den Verkaufspreis? Exportieren Sie eine Regressionstabelle, die beide Modelle beinhaltet.

25.2 Webscraping / Tidying

  1. Betrachten Sie den Wikipedia-Eintrag zum Klima in Bern: https://de.wikipedia.org/wiki/Bern#Klima. Lesen Sie die Tabelle “Monatliche Durchschnittstemperaturen und -niederschläge für Bern 1981-2010” ein. Verwenden Sie hierfür die Tools aus dem Kapitel “Datenimport aus HTML/XML”“, z.B. das Package rvest.

  2. Konzentrieren Sie sich auf die ersten drei Zeilen (Monat, Max. Temperatur, Min. Temperatur) und säubern Sie die Daten (vgl. Kapitel “Arbeiten mit Strings”), um auf folgende (oder hübschere) Tabelle zu kommen:

kable(tabelle)
Monat Max Min
Januar 2.8 -3.6
Februar 4.7 -3.1
Mrz 9.5 0.2
April 13.4 3.0
Mai 18.2 7.4
Juni 21.6 10.5
Juli 24.3 12.5
August 23.7 12.3
September 19.1 8.9
Oktober 13.8 5.4
November 7.3 0.4
Dezember 3.5 -2.3

Achten Sie auf gut lesbaren Code! Der Einreichungstermin ist der 5. Januar 2024.

Viel Erfolg!

26 Troubleshooting

26.1 Allgemeines Vorgehen bei einem Error:

  • Die meisten Fehlermeldungen haben eigentlich damit zu tun, dass irgendetwas nicht da ist: Library, Funktion, Objekt, Pfad. D.h. etwas ist nicht installiert, nicht geladen, nicht erzeugt worden, oder falsch geschrieben. Prüfen Sie, ob es Tippfehler gibt und alles geladen/erzeugt ist.
  • Fragen Sie ihren Nachbarn, ob er einen augenscheinlichen Fehler sieht. Oft ist man blind für die eigenen Tippfehler.
  • Googlen Sie den allgemeinen Teil der Fehlermeldung. I.d.R. hatte jemand vor Ihnen dasselbe Problem. Google Treffer führen häufig zu stackoverflow.com oder anderen Programmier-Hilfeforen. Dort macht es meistens Sinn:
    • Die Frage/Situation die ins Forum gepostet wurde kurz mit der eigenen Frage/Situation abzugleichen und wenn das Problem ähnlich gelagert scheint,
    • die bestbewertete Antwort zu lesen
  • Immer noch keine Lösung: Mich fragen. Am besten den Code, mit dem ich den Fehler reproduzieren kann, beilegen, sowie den Wortlaut der Fehlermeldung und Verweis auf die Zeile, die den Fehler erzeugt.

26.2 Häufigste Fehlermeldungen

  1. could not find function

Hat meistens eine dieser zwei Ursachen:

  • das Paket, das die Funktion enthält wurde nicht geladen. mit ??funktionsname kann man suchen, in welchem Paket die Funktion steckt
  • Falsch geschrieben
  1. Error in library
  • entweder falsch geschrieben oder
  • nicht installiert: install.packages("paketname"). Anführungszeichen nicht vergessen!!
  1. Error: object ‘…’ not found
  • falsch geschrieben oder
  • objekt wurde noch nicht erzeugt (Codezeile vergessen auszuführen; R wurde neu gestartet und das Objekt nicht neu erstellt).
  • überlegen, ob das Objekt direkt im Workspace liegt, oder ob es ein Objekt innerhalb eines Datensatzes ist: income -> “not found”. daten$income -> funktioniert
  1. cannot open the connection
  • Pfad/Dateiname oder Ähnliches falsch geschrieben
  • Anführungszeichen um die Pfadangabe vergessen
  • R sucht im falschen Pfad. Prüfen mit getwd(). Beim Eingeben des Pfads hilft das Drücken von TAB bei der Autovervollständigung des Pfads, was Fehler vermeidet
  • Falls der Pfad nicht bekannt ist: Den Pfad einer Datei kann man i.d.R. heruasfinden indem man die Eigenschaften der Datei (Windows: Rechtsklick; Apple: Apfel+Click) betrachtet (schauen nach “Ort” oder “Pfad” oder Ähnliches).
  1. … is not a …

R erwartet einen anderen Datentyp/Klasse als man ihm gibt.

  • prüfen ob der Vektor vom richtigen Typ ist (werden Zahlen und Buchstaben zusammengerechnet? Ist ein Vektor z.B. ein Faktor und der andere ein String?
  1. items to replace is not a multiple of / replacement has …
  • I.d.R.: Die Länge von zwei Vektoren passt nicht zusammen.
  • Prüfen, ob das was man zuweist die erwartete Länge hat, z.B. mit length() oder nrow().
  1. Error in if
  • im if-Statement steht vermutlich etwas das nicht vom Typ “logical” (d.h. TRUE oder FALSE) ist.
  • Prüfen, was genau für Daten im if-Statement verwendet werden.
  • ggf. sind Missings (NA, NaN) innerhalb der Daten.
  1. Es kommt keine Meldung aber es passiert einfach nichts

R ist nicht bereit für neue Eingaben!

  • Wenn beispielsweise ein Anführungszeichen oder eine Klammer geöffnet wurde aber nie geschlossen, dann wartet R, bis das schliessende Symbol (Klammer, Anführungszeichen) kommt. Selbst wenn man Enter drückt wird der Befehl dann nicht abgebrochen. Das erkennt man daran, dass die Konsole ein “+” zeigt statt ein “>”, d.h. das vorherige Kommandos um weitere Eingaben ergänzen will, bis das fehlende Zeichen kommt.

  • Evtl “hängt” die Session. Eine komplizierte Rechnung kann auch mal eine Weile dauern. Erkennbar i.d.R. daran, dass in der Konsole das Zeichen “>” nicht da ist sondern nur ein blinkender Coursor. Zudem erscheint rechts oben in der Konsole ein rotes Stopschild.

  • überlegen, ob R etwas sinnvolles rechnet: fertig rechnen lassen oder abbrechen mit Klick aufs Stopschild oder Escape-Taste
  • evtl. haben Sie etwas programmiert, was niemals endet (z.B. eine Endlosschleife) oder anderweitig die Ressourcen des Arbeitsgeräts in die Knie zwingt. In dem Fall abbrechen mit Escape und den Code überprüfen.