-------------------------
#  Regression-Examples  #
-------------------------

x = seq(from=-5,to=5,by=0.5)
n = length(x)
n
noise = rnorm(n)
y = 3 - 0.5 * x + noise           # Gleichung (1)
plot(x,y)

# Regression: Nehmen wir an, die Daten x und y sind gegeben, 
# wir kennen aber nicht den genauen Zusammenhang zwischen 
# x und y. Man schaut sich den plot an und vermutet einen 
# linearen Zusammenhang:
#
#   y = beta0 + beta1 * x           Gleichung (2)
#
# Die lm()-Funktion bestimmt diejenigen beta0 und beta1, 
# die am besten, sum_i (y_i^observed - y_i^model)^2 -> min, 
# zu den gegebenen Daten passt. Dabei ist hier y^observed 
# durch (1) gegeben und y^model durch Gleichung (2).

lm( y ~ x )

# der Rueckgabetyp der lm()-Funktion ist eine Liste mit 
# 12 Elementen, die default-maessig nicht alle angezeigt 
# werden.

res = lm( y ~ x )
res

summary(res)

names(res)

# structural information:
mode(res)
str(res)
class(res)

# list elements:
res[[1]]       # Listen-Elemente werden mit doppelten 
               # eckigen Klammern aufgerufen
res[1]         # ok, scheint hier auch zu gehen..
names(res)
res$coefficients
res$coef
res$c          # es gibt noch ein res$call
res$co
res$ca
res$call

# Regression-Fit:
res$fit
plot(x,y)
lines(x,res$fit,col="red")


----------------
#  Beispiel 2  #
----------------
rm(list=ls())                 # loescht alles

x = seq(from=-5,to=5,by=0.5)
n = length(x)
n
noise = rnorm(n)
y = 3 - 0.5 * x + 0.25 * x^2 + noise         
plot(x,y)

# durch Betrachten des plots vermutet man 
# einen quadratischen Zusammenhang:
#
# y = beta0 + beta1 * x + beta2 * x^2         (3)
#
# Das Modell (3) ist nicht mehr linear in x, 
# aber immer noch linear in den gesuchten 
# Modell-Parametern beta0, beta1, beta2 -> 
# wir koennen die betas mit linearer Regression, 
# mit der lm()-Funktion, bestimmen. Die 
# Syntax fuer den x^2-Term ist ein bischen 
# gewoehnungsbeduerftig:

lm( y ~ x + I(x^2) )

lm( y ~ x + x^2 )     # macht was anderes..

res = lm( y ~ x + I(x^2) )

plot(x,y)
lines(x,res$fit,col="red")  


# die lm()-Funktion nimmt automatisch immer einen
# konstanten Term mit, sie produziert automatisch 
# immer ein beta0 (auch "intercept" genannt); moechte 
# man das unterdruecken, also ein Modell 
#
# y = beta1 * x + beta2 * x^2                 (4)
#
# fitten, dann muss man folgende Syntax benutzen:

lm( y ~ -1 + x + I(x^2) )

# In diesem Fall kommt das natuerlich nicht so gut 
# hin:

res2 = lm( y ~ -1 + x + I(x^2) )
plot(x,y)
lines(x,res2$fit,col="blue")  

plot(x,y,ylim=c(0,15))
lines(x,res$fit,col="red")
lines(x,res2$fit,col="blue")