## R Ue4T2.1 dats_01_22 Konsumfunktion, Parametertests (t-, F-Test)
## Bezeichnungen siehe dat_01_22.docx 
## cpn, cpr ... Privater Konsum nominell, real (Mio Euro)
## yhn, yhr ... Haushaltseinkommen nom, real (Mio Euro)
## wln      ... Geldvermoegen nom (Mio Euro ?)
## pcd      ... Konsumdeflator 2015=100

#setwd("C:/MH/WU/LV/OEKONOMETRIE_BA/Oe1_WS23/Chp4/EXERCISES/")
setwd("C:/Users/hoersaal/Downloads/")
source("BasicStatistics_R.txt")  # Liest bereits geschriebene R-Scipts ein

# Daten einlesen
dat <- read.table("dats_01_22.csv", sep=";", dec=",", header=TRUE, na.strings = "#NV", fill = TRUE, comment.char="")

## Daten anschauen
head(dat)
tail(dat)
dim(dat)    # Zeitreihen n=47 k=9 (Variablen)
names(dat)  # Namen der Variablen

cpr <- ts(dat$cpr, start=1976)
yhr <- ts(dat$yhr, start=1976)
wlr <- ts(dat$wln/dat$pcd, start=1976)     # reales Geldvermoegen in 100 Mio 2015-Euro
wlr1 <- lag(wlr,-1)
n_cpr <- length(cpr)
year <- ts(c(1976:(1976 + n_cpr - 1)) )

## PLOT 
# 2 Reihen in einem plot mit unterschiedlichen Skalen (nicht ganz einfach)
ymax <- max(c(yhr,cpr), na.rm=T);   ymin <- min(c(yhr,cpr), na.rm=T)

par(mar=c(5, 4, 4, 5) +0.1) # Groesse des Plots, 
   # +0.1 ist der zusaetzliche Platz rechts fuer "ur"
plot.ts(cpr, ylab="yhr, cpr", col="red", axes=FALSE, ylim=c(ymin,ymax), xlim=c(1975,2025) )
lines(yhr, col="black")
axis(2,at= 50000*(2:4))
par(new=TRUE)   # naechster Plot wird daruebergelegt
plot.ts(wlr1, col="blue", axes=FALSE, ylab="")
axis(4)
mtext("wlr",side=4,col="blue",line=2) 
axis(1)
title("yhr, cpr und wlr(-1)")
legend("bottomright", legend=c("yhr","cpr","wlr(-1)"), text.col=c("black","red","blue")) 
box()


## Modell: cpr = a + b yhr + c wlr(-1) + d cpr(-1) + u 
# Datenaufbereitung
n_cpr <- length(dat$cpr)
cpr1  <- c(NA, cpr[1:(n_cpr-1)])     # muessen den Lag selbst erzeugen
wlr1  <- c(NA, wlr[1:(n_cpr-1)])     # muessen den Lag selbst erzeugen
df1   <- data.frame(year,cpr, yhr, wlr1, cpr1) 
head(df1); tail(df1)
df1   <- na.omit(df1)                # entfernen v Beobachtungen mit NA
# Schaetzung
mod <- lm(cpr ~ yhr + wlr1 + cpr1 + 1, data=df1)  
summary(mod)    


## Wo passt das Modell, wo nicht?       
# Plot v observed fitted, residuals
plot_ser_fit_res.ts(df1$cpr,mod,"cpr", 1977) 
                                
# Daten (actual), fitted, residuals
df2 <- cbind(df1,mod$fitted, mod$residuals)
tail(df2,10)