BaseURL = "http://www2.jura.uni-hamburg.de/instkrim/kriminologie/Mitarbeiter/Enzmann/Lehre/StatIKrim/"

# ------------------------------------------------------------------------------
# Aufgabe: Demonstration des Regressionseffekts, vgl. Nachtigall & Suhl (2002)
#         http://www.metheval.uni-jena.de/materialien/reports/report_2002_02.pdf
# Projekt: Statistik 1 für Kriminologen
# Autor:   Dirk Enzmann \ 12-Jan-2012

# #0
# program setup

library(foreign)
source(paste(BaseURL,"Uebungen/plot.kdnc.r",sep=""))

# ==============================================================================
# #1
# Beispieldatensatz laden:

Daten = read.dta(paste(BaseURL,"regreff.dta",sep=""))
attach(Daten)

# ------------------------------------------------------------------------------
# #2
# Grafiken und deskriptive Statistiken der Körpergröße der Mütter und Töchter:

windows(width=9.5,record=T)

mean = round(mean(groemu))
sd = round(sd(groemu),2)
N = length(groemu)
cat("Mütter: Mean = ",mean,", sd = ",sd,", N = ",N,"\n",sep="")

plot.kdnc(groemu,xlab="Körpergröße Mütter",
          title="Dichteverteilung mit Normalverteilungskurve")
dens = density(groemu)
rangey = range(dens$y)
rangex = range(dens$x)
posx = rangex[2]-(rangex[2]-rangex[1])/5
posy = rangey[2]-(rangey[2]-rangey[1])/20
liney = (rangey[2]-rangey[1])/20
text(x=posx,y=posy,paste("mean =",mean),adj=c(0,0))
text(x=posx,y=posy-liney,paste("     sd =",sd),adj=c(0,0))
text(x=posx,y=posy-2*liney,paste("       N =",N),adj=c(0,0))

mean = round(mean(groesse))
sd = round(sd(groesse),2)
N = length(groesse)
cat("Töchter: Mean = ",mean,", sd = ",sd,", N = ",N,"\n",sep="")

plot.kdnc(groesse,xlab="Körpergröße Töchter",
          title="Dichteverteilung mit Normalverteilungskurve")
dens = density(groesse)
rangey = range(dens$y)
rangex = range(dens$x)
posx = rangex[2]-(rangex[2]-rangex[1])/5
posy = rangey[2]-(rangey[2]-rangey[1])/20
liney = (rangey[2]-rangey[1])/20
text(x=posx,y=posy,paste("mean =",mean),adj=c(0,0))
text(x=posx,y=posy-liney,paste("     sd =",sd),adj=c(0,0))
text(x=posx,y=posy-2*liney,paste("       N =",N),adj=c(0,0))

# ------------------------------------------------------------------------------
# #3
# Demonstration der Regressionseffekte bei Auswahl von Extremgruppen:

groemu_z = scale(groemu)
groeto_z = scale(groesse)

# Kleine Töchter großer Mütter bzw. große Töchter kleiner Mütter:
meanmumug = round(mean(groemu_z[groemu_z > 1]),2)
sdmumug = round(sd(groemu_z[groemu_z > 1]),3)
Nmumug = length(groemu_z[groemu_z > 1])
meantomug = round(mean(groeto_z[groemu_z > 1]),2)
sdtomug = round(sd(groeto_z[groemu_z > 1]),3)
Ntomug = length(groeto_z[groemu_z > 1])

meanmumuk = round(mean(groemu_z[groemu_z < -1]),2)
sdmumuk = round(sd(groemu_z[groemu_z < -1]),3)
Nmumuk = length(groemu_z[groemu_z < -1])
meantomuk = round(mean(groeto_z[groemu_z < -1]),2)
sdtomuk = round(sd(groeto_z[groemu_z < -1]),3)
Ntomuk = length(groeto_z[groemu_z < -1])

cat("Größe Mütter  (z-Mütter >  1): Mean = ",meanmumug,", sd = ",sdmumug,", N = ",Nmumug,"\n",sep="")
cat("Größe Töchter (z-Mütter >  1): Mean = ",meantomug,", sd = ",sdtomug,", N = ",Ntomug,"\n",sep="")
cat("Größe Mütter  (z-Mütter < -1): Mean = ",meanmumuk,", sd = ",sdmumuk,", N = ",Nmumuk,"\n",sep="")
cat("Größe Töchter (z-Mütter < -1): Mean = ",meantomuk,", sd = ",sdtomuk,", N = ",Ntomuk,"\n",sep="")

gewto_z = scale(gewicht)
# Geringes Gewicht großer Töchter und umgekehrt:
meangewg = round(mean(gewicht[groeto_z > 1],na.rm=T),2)
sdgewg = round(sd(gewicht[groeto_z > 1],na.rm=T),3)
Ngewg = length(gewicht[groeto_z > 1 & is.na(gewicht)==F])
meangewzg = round(mean(gewto_z[groeto_z > 1],na.rm=T),2)
sdgewzg = round(sd(gewto_z[groeto_z > 1],na.rm=T),3)
Ngewzg = length(gewto_z[groeto_z > 1 & is.na(gewto_z)==F])

meangewk = round(mean(gewicht[groeto_z < -1]),2)
sdgewk = round(sd(gewicht[groeto_z < -1]),3)
Ngewk = length(gewicht[groeto_z < -1])
meangewzk = round(mean(gewto_z[groeto_z < -1]),2)
sdgewzk = round(sd(gewto_z[groeto_z < -1]),3)
Ngewzk = length(gewto_z[groeto_z < -1])

cat("  Gewicht  (z Größe >  1): Mean = ",meangewg, ", sd = ",sdgewg, ", N = ",Ngewg, "\n",sep="")
cat("z(Gewicht) (z Größe >  1): Mean = ",meangewzg,", sd = ",sdgewzg,", N = ",Ngewzg,"\n",sep="")
cat("  Gewicht  (z Größe < -1): Mean = ",meangewk, ", sd = ",sdgewk, ", N = ",Ngewk, "\n",sep="")
cat("z(Gewicht) (z Größe < -1): Mean = ",meangewzk,", sd = ",sdgewzk,", N = ",Ngewzk,"\n",sep="")

# Kleine Mütter großer Töchter und große Mütter kleiner Töchter:
meanmutog = round(mean(groemu_z[groeto_z > 1]),2)
sdmutog = round(sd(groemu_z[groeto_z > 1]),3)
Nmutog = length(groemu_z[groeto_z > 1])
meantotog = round(mean(groeto_z[groeto_z > 1]),2)
sdtotog = round(sd(groeto_z[groeto_z > 1]),3)
Ntotog = length(groeto_z[groeto_z > 1])

meanmutok = round(mean(groemu_z[groeto_z < -1]),2)
sdmutok = round(sd(groemu_z[groeto_z < -1]),3)
Nmutok = length(groemu_z[groeto_z < -1])
meantotok = round(mean(groeto_z[groeto_z < -1]),2)
sdtotok = round(sd(groeto_z[groeto_z < -1]),3)
Ntotok = length(groeto_z[groeto_z < -1])

cat("Größe Mütter  (z-Töchter >  1): Mean = ",meanmutog,", sd = ",sdmutog,", N = ",Nmutog,"\n",sep="")
cat("Größe Töchter (z-Töchter >  1): Mean = ",meantotog,", sd = ",sdtotog,", N = ",Ntotog,"\n",sep="")
cat("Größe Mütter  (z-Töchter < -1): Mean = ",meanmutok,", sd = ",sdmutok,", N = ",Nmutok,"\n",sep="")
cat("Größe Töchter (z-Töchter < -1): Mean = ",meantotok,", sd = ",sdtotok,", N = ",Ntotok,"\n",sep="")

# ==============================================================================