### Eigene Funktion zur Durchfuehrung des Hausman-Test
# Argumente:
# - model_re: Random effects Modell (mittels lme geschaetzt)
# - model_fe: Fixed effects Modell (mittels lm geschaetzt)
# - fixed_effects: Vektor mit den Namen der fixed effects aus dem Modell,
#                  bzgl. welcher getestet werden soll
#                  Bsp.: c("x1","x2","x3")
# - alpha: Signifikanzniveau (default 0.05)
hausman.test <- function(model_re, model_fe, fixed_effects, alpha = 0.05) {
  
  # Komponenten fuer Test extrahieren
  beta_re <- summary(model_re)$tTable[,"Value"]
  beta_fe <- summary(model_fe)$coefficients[,"Estimate"]
  var_re <- vcov(model_re)
  var_fe <- vcov(model_fe)
  
  # Komponenten auf die angegebenen fixed_effects einschraenken
  # (Vektoren (und Zeilen/Spalten der Matrizen) werden gleichzeitig in gleiche Reihenfolge sortiert)
  beta_re <- beta_re[fixed_effects]
  beta_fe <- beta_fe[fixed_effects]
  var_re <- var_re[fixed_effects, fixed_effects]
  var_fe <- var_fe[fixed_effects, fixed_effects]
  
  
  # Test durchfuehren
  T <- t(beta_fe - beta_re) %*% solve(var_fe - var_re) %*% (beta_fe - beta_re)
  p <- length(beta_fe)
  chisq_alpha <- qchisq(1 - alpha, df = p) # (1-alpha)-Quantil der Chi^2-Verteilung mit df = p
  testentscheidung <- T > chisq_alpha
  print(paste0("Testentscheidung: ", ifelse(testentscheidung, "H_1", "H_0"),
               " (p = ", round(1 - pchisq(T, df = p), 4), ")"))
}