---
title: "Multiple Regression"
author: "Übungsaufgabe zu Design, Analyse, Dokumentation SoSe 2023"
bibliography: referenzen.bib
lang: en
format: 
  pdf:
    include-in-header:
      text: |
       \usepackage[font=small,format=plain,labelfont=bf, labelsep=period,justification=justified,singlelinecheck=false]{caption}
---

Grundlage und Inspiration für diese Übungsaufgabe ist die Studie von @covi1976. 
Ziel ist es, mithilfe von partiellen Korrelationsanalysen und einer multiplen Regressionsanalyse
zu quantifizieren, welche Einfluss das Patient:innenalter und die Qualität der 
Therapeut:innen-Patient:innenbeziehung auf den Erfolg einer Psychotherapie bei 
Depression haben. Dazu betrachten wir an @covi1976  beispielgebend nur die Gruppe 
der Gruppentherapieintervention und fokussieren auf simulierte Pre-Post *Profiles of Mood (POMS)* 
Depressionsselbsteinschätzungsdifferenzwerte als fiktives Ergebnismaß der Studie. 

## Datensatz {-}

Der Datensatz `7-Multiple-Regression.csv` enthält als erste Spalte simulierte 
Pre-Post POMS Depressionsselbsteinschätzungsdifferenzwerte (`dPOMS`) für insgesamt 
$n = 100$ Patient:innen, als zweite Spalte das Patient:innenalter (`AGE`) und 
als dritte Spalte die Qualität der Patient:innen-Therapeut:innenbeziehung 
basierend auf dem Barrett-Lennard Relationship Inventory (`BRI`). Alle Variablendaten 
liegen dabei in *normalisierter* Form, d.h. mit Stichprobenmittel 0 und 
Stichprobenvarianz 1 vor. @tbl-bdi zeigt exemplarisch die Daten der ersten zehn Patient:innen.

```{r, echo = F}
# Modellformulierung
set.seed(1)                                                                     # Ergebnisreproduzierbarkeit
library(MASS)                                                                   # Multivariate Normalverteilung
n         = 100                                                                 # Anzahl Patient:innen
X         = round(mvrnorm(n, c(0,0),diag(2)), digits = 1)                       # Kovariatenrealisierungen 
beta      = matrix(c(1,1), nrow = 2)                                            # wahrer, aber unbekannter, Betaparameter
sigsqr    = 1                                                                   # wahrer, aber unbekannter, Varianzparameter
y         = mvrnorm(1, X%*%beta, sigsqr*diag(n))                                # Datenrealisierung
D         = data.frame(dPOMS = y,                                               # dPOMS Werte
                       AGE   = X[,1],                                           # AGE Regressorwerte
                       BRI   = X[,2])                                           # BRI Regressorwerte
fname   = file.path("7-Multiple-Regression.csv")                                # Dateiname
write.csv(D, file = fname, row.names = FALSE)                                   # Speicherung 
```

\footnotesize
```{r echo = F, eval = T}
#| label: tbl-bdi
#| tbl-cap : "Daten der ersten 10 Patient:innen des Datensatzes"
library(knitr)                                                                  # knitr für Tabellen
fname       = "7-Multiple-Regression.csv"                                       # Dateiname
D           = read.table(file.path(fname), sep = ",", header = TRUE)            # Laden des Datensatzes
kable(D[1:10,], digits = 2, align = "c")                                        # erste 10 Patient:innen
``` 

\newpage
\normalsize
## Programmieraufgaben {-}
\noindent 1. Bestimmen Sie die paarweisen Korrelationen der Variablen 
`dPOMS`, `AGE` und `BRI`. Bestimmen Sie weiterhin die partiellen Korrelationen 
von `dPOMS` und `AGE` gegeben `BRI` und von `dPOMS` und `BRI` gegeben `AGE` mithilfe
der **R** Funktion `ppcor()`. Bestimmen Sie schließlich den Betaparameterschätzer 
in einer multiplen Regressionsanalyse mit abhängiger Variable `dPOMS` und 
unabhängigen Variablen `AGE` und `BRI`, einmal mithilfe der ALM Standardformel
und einmal auf Grundlage der partiellen Korrelationen.

\footnotesize
```{r, warning = F, echo = F}
fname      = file.path(getwd(), "7-Multiple-Regression.csv")                    # Dateipfad     
D          = read.table(fname, sep = ",", header = TRUE)                        # Datensatzeinlesen
y          = D$dPOMS                                                            # Abhängige Variable
n          = length(y)                                                          # Anzahl Datenpunkte
X          = matrix(c(rep(1,n), D$AGE, D$BRI), nrow = n)                        # Designmatrix
p          = ncol(X)                                                            # Anzahl Parameter
beta_hat   = solve(t(X) %*% X) %*% t(X) %*% y                                   # Betaparameterschätzer
eps_hat    = y - X %*% beta_hat                                                 # Residuenvektor
sigsqr_hat = (t(eps_hat) %*% eps_hat) /(n-p)                                    # Varianzparameterschätzer

# Betaparameterschätzer aus partiellen Korrelationen und Korrelationen
library(ppcor)                                                                  # partielle Korrelationentoolbox
y12        = cbind(y,X[,2:3])                                                   # y,x_1,x_2 Matrix
bars       = apply(y12, 2, mean)                                                # Stichprobenmittel
s          = apply(y12, 2, sd)                                                  # Stichprobenstandardabweichungen
r          = cor(y12)                                                           # Stichprobenkorrelationen
pr         = pcor(y12)                                                          # partielle Stichprobenkorrelationen
pr         = pr$estimate                                                        # partielle Stichprobenkorrelationen
beta_hat_1 = pr[1,2]*sqrt((1-r[1,3]^2)/(1-r[2,3]^2))*(s[1]/s[2])                # \hat{\beta}_1
beta_hat_2 = pr[1,3]*sqrt((1-r[1,2]^2)/(1-r[3,2]^2))*(s[1]/s[3])                # \hat{\beta}_2
beta_hat_0 = bars[1] - beta_hat_1*bars[2] - beta_hat_2*bars[3]                  # \hat{\beta}_0

# Ausgabe
cat("Korrelationen  r(dPOMS,AGE),r(dPOMS,BRI),r(AGE,BRI)        :"  , round(c(r[1,2],r[1,3],r[2,3]), digits = 3),
     "\nPartielle Korrelationen r(dPOMS,AGE|BRI), r(dPOMS,BRI|AGE) :", round(c(pr[1,2],pr[1,3]), digits = 3),
     "\nbeta_hat ALM Schätzer                                      :", round(beta_hat, digits = 3),
     "\nbeta_hat aus partieller Korrelation                        :", round(c(beta_hat_0,beta_hat_1,beta_hat_2), digits = 3))
```

\normalsize

\noindent 2. Visualisieren den Datensatz zusammen mit der durch den Betaparameterschätzer
definierten Regressionsebene. Die Abbildung sollte in etwa aussehen wie @fig-mrg-fit.

```{r, eval = F, echo = F, warning=FALSE, message=FALSE}
# Dateneinlesen
fname      = file.path(getwd(), "7-Multiple-Regression.csv")
D          = read.table(fname, sep = ",", header = TRUE)        
library(scatterplot3d)
graphics.off()
pdf(
file        = file.path("7-Multiple-Regression-Abbildung.pdf"),
width       = 5,
height      = 5)
spd = scatterplot3d(
D$AGE,
D$BRI,
D$dPOMS,
pch         = 21,
color       = "gray",
bg          = "black",
type        = "h",
xlab        = "AGE",
ylab        = "BRI",
zlab        = "dPOMS",
angle       = 75)
alm         = lm(dPOMS ~ AGE + BRI, data = D)
spd$plane3d(alm)
dev.off()
``` 

![Visualisierung des Datensatzes und Regressionsebene](7-Multiple-Regression-Abbildung){#fig-mrg-fit fig-align="center" width=70%}  


\newpage
## Dokumentation

Bitte beachten Sie bei der Erstellung Ihre Dokumentation folgende Vorgaben und 
orientieren Sie sich in der Darstellung Ihrer datenanalytischer Ergebnisse an
den Empfehlungen des [APA Publication Manuals 7th Edition](https://apastyle.apa.org/products/publication-manual-7th-edition), insbesondere Kapitel 6. 

### Einleitung {-}

Stellen Sie kurz die Ausgangsfrage der Studie von @covi1976 und speziell der vorliegenden 
Übungsaufgabe dar und erläutern Sie kurz die Profiles of Mood State (POMS) und das Barrett-Lennard Relationship Inventory (BRI) als Maße für Depressivität und die Qualität von Patient:innen-Therapeut:innenbeziehungen.

### Methoden {-}

Erläutern Sie kurz die Begriffe der partiellen Korrelation und der multiplen Regression.


### Resultate {-}

Reportieren Sie die von Ihnen in Programmieraufgabe 1 bestimmten Statistiken. Diskutieren
Sie dabei die wechselseitigen Abhängigkeiten der drei betrachteten Variablen und 
gehen Sie dabei insbesondere auf die Korrelation der `AGE` und `BRI` Variablen und ihre
jeweiligen partiellen Korrelationen mit der `dPOMS` Variable ein. Erläutern Sie
die Lage und Orientierung der Regressionsebene vor dem Hintergrund der in 
Programmieraufgabe 1 bestimmten Werte des ALM Betaparameterschätzers.

### Schlußfolgerung

Fassen Sie die von Ihnen erstellte Dokumentation in drei Sätzen zusammen.

## Referenzen 
\footnotesize