#################################################
### Vorbereitungen / Wiederholung
#################################################

### Datensatz laden und überprüfen (vorher bei Bedarf das notwendige Paket AER installieren)

library(AER)
data(CPS1985)
attach(CPS1985)


### Was ist bisher bekannt?

# High-Level-Grafikfunktionen
pairs(CPS1985[, 1:4])



# Konfigurierbarkeit von High-Level-Grafikfunktionen direkt
?pairs
pairs(CPS1985[, 1:4], labels = c("Lohn", "Ausbildung", "Berufserfahrung", "Alter"))



# Konfigurierbarkeit von High-Level-Grafikfunktionen über par()
par(bg = "red")
pairs(CPS1985[, 1:4], labels = c("Lohn", "Ausbildung", "Berufserfahrung", "Alter"))

par(bg = "transparent")		# Ursprungszustand wiederherstellen
par()					# Ursprungszustand einer globalen Grafikoption abfragen



# Grafik-Devices
par(mfrow = c(1, 2))	# 4 Graphen pro Device, in 2 Zeilen und 2 Spalten
  plot(wage ~ education)
  plot(log(wage) ~ education)
par(mfrow = c(1, 1))





#################################################
### Low-Level-Grafikfunktionen
#################################################

### Plot zusammensetzen/ergänzen

plot.default(x = c(0, 5), y = c(0,3))	# leeren Plot erstellen: sehr aufwändig
							# deshalb meist:
							# vorhandenen Graphen modifizieren
plot(education, wage, pch = 20)		# nun kann man Low-Level-Elemente hinzufügen, siehe:
axis(side = 4)					# Achse auf rechter Seite (hier nicht nötig)
rug(wage, side = 2)				# Striche auf linker Seite, die einen Eindruck der Verteilung vermitteln
							# Achsennummerierung der Graphen im Uhrzeigersinn (mit side = 1 als Abszisse)

### Geraden zeichnen

reg1 <- lm(wage ~ education)			# Regressionsgerade berechnen
abline(reg1, lwd = 2, col = "green")	# lwd = 1 wäre die normale Dicke, lwd = 2 macht die Objekte doppelt so dick.
abline(reg1, lwd = 2, col = 2)		# Beachte: col akzeptiert auch Zahlen der Farben (1 ist schwarz, 2 ist rot, usw.)
							# Vergleiche: Section "Color Specification" in ?par

# Die Funktion abline benötigt als Argumente entweder Achsenabschnitt und Steigung oder für horizontale (vertikale)
# Linien die Argumente h= (bzw. v=), vergleiche:
abline(h = mean(wage), v = mean(education))	# Zur Erinnerung: die Regressionsgerade verläuft immer durch die Mittelwerte von x und y.



### Strecken zeichnen

# Residuum für den Punkt education = 5, wage = 14 einzeichnen.
lines(x = c(5, 5), y = c(14, fitted(reg1)[education == 5]), col = "darkgray", lwd = 3)
# Beachte: x = c(x-Wert Startpunkt, x-Wert Endpunkt), y analog
# Man kann auch Linien aus mehr als zwei Punkten einzeichnen, siehe:
lines(x = c(2, 7, 10), y = c(20, 30, 25))



### Flächen zeichnen

# Residuenquadrat für Beobachtung mit education = 3
polygon(x = c(3, 3, 3 + resid(reg1)[education == 3], 3 + resid(reg1)[education == 3]), 
	y = c(wage[education == 3], fitted(reg1)[education == 3], fitted(reg1)[education == 3], wage[education == 3]))



### Punkte zeichnen

plot(education[gender == "male"], wage[gender == "male"], pch = 20, col = "blue")		# neuer Plot (Männer)
points(education[gender == "female"], wage[gender == "female"], pch = 20, col = "red")	# Frauen
summary(wage)	# Beachte: Die Ordinate geht nur bis etwa 25, aber der maximale Lohn liegt bei 44.5, deshalb:

plot(education[gender == "male"], wage[gender == "male"], pch = 20, col = "blue", 
	xlim = range(education), ylim = range(wage),
	xlab = "Ausbildung", ylab = "Lohn")		# Beschriftung anpassen
points(education[gender == "female"], wage[gender == "female"], pch = 20, col = "red")




### Weitere Informationen hinzufügen

# Regressionsgeraden für Männer und Frauen
reg1.male	<- update(reg1, data = CPS1985[gender == "male", ])
reg1.female	<- update(reg1, data = CPS1985[gender == "female", ])
abline(reg1.male, col = "blue", lwd = 2)
abline(reg1.female, col = "red", lwd = 2)



# Rugs
rug(wage[gender == "male"], side = 2, col = "blue")
rug(wage[gender == "female"], side = 4, col = "red")



### Legende

legend("topleft", 
	c("Regressionsgerade", "Männer", "Frauen"), 
	lty = c(1, NA, NA), 
	pch = c(NA, 20, 20), 
	lwd = c(2, 1, 1), 
	col = c("black", "blue", "red"))



### Pfeile einzeichnen

arrows(x0 = 10, y0 = 40, x1 = 13.8, y1 = 44.2, lwd = 2, col = "green")



### Text hinzufügen

# als Überschrift
title("Lohn versus Ausbildung")

# im Diagramm
text(x = 11.5, y = 39, labels = paste("Beobachtung Nr. ", which(wage > 40), sep = ""), cex = 1.5)



### Mathematische Beschriftungen (kurzes Simulationsbeispiel)

?plotmath	# Erläuterung der mathematischen Beschriftung

set.seed(42)

# DGP definieren
n		<- 100				# Beobachtungen
x		<- rnorm(n, mean = 0, sd = 1)	# Regressor
beta_0	<- 2					# Konstante
beta_1	<- 1					# Steigungsparameter

R 		<- 200				# Wiederholungen

# Matrix zum speichern der geschätzten Parameter
beta_dach	<- array(NA, dim = c(R, 2), dimnames = list(NULL, c("beta_0_dach", "beta_1_dach")))

# Schleife
for(r in 1:R){
  u			<- rnorm(n, mean = 0, sd = 1)	# Fehlerterm
  y			<- beta_0 + beta_1*x + u	# DGP
  reg			<- lm(y ~ x)			# Modell schätzen
  beta_dach[r, ]	<- reg$coefficients		# Geschätzte Parameter speichern
}



plot(density(beta_dach[, 2]), xlim = range(beta_dach[, 2]), lwd = 2, col = "blue", 
	xlab = expression(hat(beta)[1]), ylab = "Dichte", main = "")
lines(seq(from = min(beta_dach[, 2]), to = max(beta_dach[, 2]), length.out = 50), 
	dnorm(seq(from = min(beta_dach[, 2]), to = max(beta_dach[, 2]), length.out = 50), 
		mean = mean(beta_dach[, 2]), sd = sd(beta_dach[, 2])),
	col = "red", lwd = 2)
text(1, 1, 
	expression(f(x) == frac(1, sigma ~~ sqrt(2*pi)) ~~ e^{-frac((x - mu)^2, 2*sigma^2)}), 
	cex = 1.5, col = "red")



### Weiterhin nützlich

?segments
?grid
?mtext


#################################################
### Aufgaben I
#################################################

# 1. Lesen Sie in Ligges 2008 das Kapitel 8.1.

# 2. Versuchen Sie einen Christbaum in einem Plot mit möglichst wenigen Zeilen Code darzustellen.

# 3. Analysieren Sie einen Zusammenhang zwischen den Variablen im Datensatz CPS1985, der Ihnen wichtig erscheint 
#    mittels einer aussagekräftigen Graphik, in der Sie mindestens 4 verschiedene Low-Level-Grafikfunktionen 
#    verwenden.

# 4. Erläutern Sie die ökonometrische Intuition, die hinter dem Beispiel zur mathematischen Beschriftung steckt.
#    Verwenden Sie dazu bei Bedarf ein paar ausgewählte Grafiken.

# 5. Erzeugen Sie ein leeres Plot-Fenster, in welches Sie die beiden Formeln für die KQ-Schätzer von Achsenabschnitt
#    und Steigung im linearen Einfachregressionsmodell eintragen.