# --------- Vectorizing Computations  -------
# Programmieren ohne explizite Schleife

rm(list = ls())

# ----------------------------------------
# 6.1 Vectorization

set.seed(1234)
x <- rnorm(5)
y <- runif(5)

# Vektorisierte Berechnung (implizite Schleife)
x + y 

# statt expliziter Schleife
z <- numeric(length(x))
for ( i in seq(along.with = x) ) {
  z[i] <- x[i] + y[i]
}
z

# Vektorisierte Berechnungen sind
# - kompakter
# - besser lesbar
# - meist schneller (bei größen Eingaben)
set.seed(1234)
x <- rnorm(10^6)
y <- runif(10^6)

# Implicit loop: interne Schleife, in Maschinencode programmiert (daher schneller)
system.time(x + y)

# Explicit loop: explizite Schleife, in R programmiert
z <- numeric(length(x))
system.time(for ( i in seq(along.with = x) ) z[i] <- x[i] + y[i])

# -> wannimmer möglich, sollte in R vektorisiert programmiert werden

# Vorsicht: in R werden "Recycling rules" verwendet.
# Wenn mit zwei Vektoren unterschiedlicher Länge gearbeitet wird, 
# und die Länge des einen Vektors ein Vielfaches des anderen ist,
# dann wird der kürzere Vektor wiederholt.
c(1, 2, 3) + 1:6
# Nur wenn die Länge des einen Vektors kein Vielfaches des anderen ist,
# wird eine Fehlermeldung ausgegeben.
c(1, 2, 3) + 1:5

# ----------------------------------------
# 6.2 Implicit loops
# Beispiele - hier sind implizite Schleifen direkt enthalten

#--------
# Operatoren 
?Arithmetic
?Comparison

set.seed(1234)
x <- sample(100, 5)
((x %% 2) == 0)

#--------
# Statistische Maßzahlen
mean(x)    # Summe als implizite Schleife, geteilt durch n

sum(x) / length(x)

# dasselbe mit expliziter Schleife
z <- 0
for ( i in seq(along.with = x) ) {
  z <- z + x[i]
}
z <- z / length(x)
z

# auch möglich für Datensätze - Zeilen- und Spalten-Summen (und Mittelwerte)
data("cars", package = "datasets")
?cars
str(cars)
colMeans(cars)  # implizite Schleife

?rowMeans
?rowSums
?colSums

#--------
# Subsets - Extraktion bzw. Ersetzungen von Teilmengen
set.seed(1234)
x <- letters[1:10]
ind <- sample(100, 10)

# Teilmenge der Buchstaben, für die ind gerade ist
x[ind %% 2 == 0]
# bzw. mit which()
x[which(ind %% 2 == 0)]

# als explizite Schleife
z <- character()
for ( i in seq(along.with = x) ) {
  if ( ind[i] %% 2 == 0 ) {
    z <- c(z, x[i])
  }
}
z

# oder mit ifelse (aber längerer Vektor mit zusätzlichen NA-Einträgen)
ifelse(test = ind %% 2 == 0, yes = x, no = NA)

#--------
# Matrizenrechnung
# für Matrizen (zweidimensional) - matrix
# und Arrays (drei- und mehrdimensional) - array

# z.B. Matrixmultiplikation, Transponierung, Bilden von Teilmengen

# Transponieren mit impliziter Schleife
set.seed(1234)
A <- matrix(sample(10), ncol = 2)
A
t(A)

# Transponieren mit expliziter Schleife
n <- nrow(A)
m <- ncol(A)
At <- matrix(NA, nrow = m, ncol = n)
for ( i in seq(length.out = n) ) {
  for ( j in seq(length.out = m) ) {
    At[j, i] <- A[i, j]
  }
}
At


# ----------------------------------------
# 6.3 The apply-family
# hiermit kann man für alle weiteren Funktionen implizite Schleifen programmieren
# wiederholte Anwendung einer Funktion für
# - alle Elemente eines Vektors
# - alle Zeilen / Spalten einer Matrix
# - alle Dimensionen eines Arrays

# Vorteile von apply gegenüber expliziter Schleife
# - kompaktere Berechnung
# - Parallelisierbar (mit Paket plyr)
# - oft schneller

# Beispiele

#------
# Funktion angewendet auf eine Liste bzw. einen Vektor
ab <- function(x) {
  if ( x < 10 ) {
    "a"
  } else {
    "b"
  }
}

# anwendbar für Skalare
ab(4)
ab(20)

# nicht korrekt anwendbar für Vektoren
ab(5:15)

# anwendbar für Vektoren mit 
?lapply
# gibt immer eine Liste zurück
lapply(X=5:15, FUN=ab)

?sapply
# versucht, das Ergebnis zu vereinfachen 
# (zu einem Vektor, einer Matrix, einem Array)
sapply(X=5:15, FUN=ab)



# Liste mit Ergebnissen von verschiedenen Durchgängen eines Experiments
set.seed(1234)
l <- list(exp1 = rnorm(100),
          exp2 = rnorm(100),
          exp3 = rnorm(100))

# Mittelwerte (von Hand berechnet)
c(mean(l$exp1),
  mean(l$exp2),
  mean(l$exp3))

# mit expliziter Schleife
ms <- numeric(length = length(l))
for ( i in seq(along.with = l) ) {
  ms[i] <- mean(l[[i]])
}
ms

# einfacher mit
sapply(l, mean)

# Unterschied zwischen lapply und sapply 
lapply(l, quantile)  # gibt Liste zurück
sapply(l, quantile)  # gibt Matrix zurück 

# multivariate Version von sapply
?mapply
mapply(rep, 1:4, times=4:1)
mapply(FUN=quantile, l, probs=list(c(0,0.5,0.75,1), c(0.2,0.4), c(0.1,0.3,0.6)))

# wiederholte Auswertung von Ausdrücken mit
?replicate
replicate(n=5, expr=rnorm(10,0,1))

#------
# Funktion angewendet auf eine Matrix, ein Array oder einen Datensatz
data("cars", package = "datasets")

?apply
# Maximum pro Zeile bzw. Spalte
apply(X=cars, MARGIN=1, FUN=which.max)
apply(X=cars, MARGIN=2, FUN=which.max)
# MARGIN=1 für Zeilen
# MARGIN=2 für Spalten
# MARGIN=c(1,2) für Zellen


# ----------------------------------------
# 6.4 Anonymous functions

# Funktion ohne Name
function(x) {
  x + 1
}

# Anwendung innerhalb von apply-Befehlen
sapply(5:15,
       function(x) {
         if ( x < 10 ) {
           "a"
         } else {
           "b"
         }
       })


# ----------------------------------------
# 6.5 apply-family gems
# Vereinfachung des Ergebnisses bzw. Zugriff auf Teile des Ergebnisses

ds <- lapply(1:10, function(i) data.frame(a = runif(4), b = gl(2, 2), c=runif(4)))
# ds <- replicate(10, data.frame(a = runif(4), b = gl(2, 2), c=runif(4)), simplify=FALSE) # alternativ 
str(ds, 1)

# Vereinfachung des Ergebnisses - in ein großes Objekt
# do.call(, lapply())

# Beispiel: Liste von Datensätzen in einen großen Datensatz umwandeln
?do.call
d <- do.call(what=rbind, args=ds)
str(d)


# Zugriff auf Teile einer Liste
# lapply(, "[[", )

# Beispiel: Spalte b für jeden Listeneintrag auswählen
bs  <- lapply(ds, "[[", "b")
# lapply(ds, "[[", c("b","c")) # funktioniert nicht 

# daher mit anonymous functions
bs2 <- lapply(1:10, function(j) ds[[j]][,c("b","c")])
bs3 <- lapply(1:10, function(j) ds[[j]][c("b","c")])


bs4 <- lapply(1:10, function(j) ds[[j]]["b"])                # Liste aus Data-frames 
bs5 <- lapply(1:10, function(j) ds[[j]][,"b"])               # Liste aus Vektoren 
bs6 <- lapply(1:10, function(j) ds[[j]][, "b", drop=FALSE])  # Liste aus Data-frames 

# auch kompliziertere Operationen/Selektionen möglich 
bs7 <- lapply(1:10, function(j) ds[[j]][2,"b"])  

# ----------------------------------------
# 6.6 Functional programming

# Funktionale Programmierung: Berechnungen als Auswertung von Funktionen

# higher-order functions - tun mindestens eines der beiden folgenden
# - nehmen eine oder mehrere Funktionen als Input
# - geben eine Funktion zurück

# -> apply-Funktionen sind higher-order functions
# -> die meisten anderen (wie mean()) sind first order functions

# weitere higher-order functions in R
?Map