#---------------------------------------------------------------------
# Thema:   Einfache lineare Regression
# Autorin: Lisa-Marie Joerißen
#
#                               ÜBUNG
#---------------------------------------------------------------------
library(stats)

#---------------------------Datensimulation---------------------------

# Seed setzen
set.seed(5)

# Simulationsparameter
n <- 21                   # Proband:innnen pro Gruppe
mu <- c(32, 29, 26)       # Erwartungswertparameter
sigsqr <- c(5, 12, 27)    # Varianzparameter

# Datensimulation als Nv mit 4 Nachkommastellen
DS <- data.frame(
  Int = round(rnorm(n, mu[1], sqrt(sigsqr[1])), 4),
  Mid = round(rnorm(n, mu[2], sqrt(sigsqr[2])), 4),
  Ter = round(rnorm(n, mu[3], sqrt(sigsqr[3])), 4)
)
View(DS)

#------------------------------Aufgabe 1:-------------------------------
# Berechnet Steigungs- und Offsetparameter der simulierten Datenwerte.
#-----------------------------------------------------------------------

# leerer Dateframe R für Results der 21 Proband:innen
R <- data.frame(
  ymean = rep(NaN, 21),
  cxy = rep(NaN, 21),
  beta_1_hat = rep(NaN, 21),
  beta_0_hat = rep(NaN, 21)
)

# Vektor x-Werte
x <- c(Int = 1, Mid = 2, Ter = 3)

#-------------------------Stichprobenstatistiken----------------------

# Stichprobenmittel der x-Werte
x_bar <- mean(x)

# Stichprobenmittel der y-Werte aller Zeilen
R$ymean <- rowMeans(DS)

# Stichprobenvarianz der x-Werte
s2x <- var(x)

# Stichprobenkovarianz der y-Werte einer Zeile und x-Werte
for (i in 1:21) {
  R$cxy[i] <- cov(as.numeric(DS[i, 1:3]), x)
}

#----------------------Ausgleichsgeradenparameter---------------------

# HINWEIS: Hier Code einfügen!

# Berechnung der Steigungsparameter




# HINWEIS: Hier Code einfügen!

# Berechnung der Offset Parameter



#---------------------Berechnung mit Funktion lm()--------------------

# Leerer Dataframe für Parameter
lm <- data.frame(
  lm_beta_1_hat = rep(NaN, 21),
  lm_beta_0_hat = rep(NaN, 21)
)

# Schleife über die Zeilen des Dataframes DS
for (i in 1:21) {
  lm_result <- lm(formula = as.numeric(DS[i, 1:3]) ~ x) # Lineare Regression
  lm$lm_beta_1_hat[i] <- coef(lm_result)[2] # Steigungsparameter extrahieren
  lm$lm_beta_0_hat[i] <- coef(lm_result)[1] # Offsetpaparameter extrahieren
}
View(lm)

#------------------------------Aufgabe 2:------------------------------
# Fragespiel: was machen die Funktionen/Argumente?
#----------------------------------------------------------------------

#------------------Visualisierungsinitialisierung----------------------

#
par(mfcol = c(2, 2),          #
  pin = c(2, 2),              #
  cex.lab = 1.25,             #
  cex.axis = 1.25,            #
  mgp = c(2.5, 1, 0),         #
  oma = c(2, 2, 2, 2),        #
  font.main = 1               #
)

#
for (i in 1:4) {
  plot(x, DS[i, 1:3],
    pch = 16,                 #
    cex = 1.25,               #
    xlab = "Therapieverlauf",
    ylab = "BDI",
    xlim = c(0.5, 3.5),       #
    ylim = c(15, 37),         #
    main = (paste("Patient:in ", i,
                  ", Steigungsparameter = ", round(lm$lm_beta_1_hat[i], 1))),
    xaxt = "n",               #
    bty = "l",                #
  )

  #
  axis(1, at = 1:3, labels = c("Int", "Mid", "Ter"))

  #
  abline(
    coef = c(lm$lm_beta_0_hat[i], lm$lm_beta_1_hat[i]),
    lty = 1,
    col = "gray"
  )

  if (i == 1) {
    legend("bottomleft", c(("Messwerte"), ("Ausgleichsgerade")),
      lty = c(0, 1),          #
      pch = c(16, NA),        #
      bty = "n",              #
      col = c("black", "gray"),
      inset = c(0.02, 0.05),  #
      yjust = 0.5,            #
      xjust = 0,              #
    )
  }
}

#---------------------------Speicherung als PDF-----------------------

# HINWEIS: Hier Pfad einfügen!

dateipfad <- "hier eigenen Pfad einfügen.pdf"
dev.copy2pdf(
  file = file.path(dateipfad),
  width = 8,
  height = 8,
)