Einleitung
Vielleicht kam Ihnen aus der Welt der Verschwörungstheorien mal zu Ohren, dass die Anzahl der COVID-Erkrankungen mit der Anzahl der 5G-Tower zusammenhängt. Aber wussten Sie, dass auch der Konsum von Eiscreme und die Anzahl der Morde in New York, oder die Anzahl von Nicolas-Cage-Filmauftritten mit der Anzahl weiblicher Redakteure beim Harvard Law Review positiv korreliert sind?
Natürlich wissen Sie bereits, dass solch eine Korrelation kein Beweis dafür ist, dass 5G-Strahlung für COVID-Erkrankungen kausal verantwortlich ist, oder dass Eiskonsum zu einer erhöhten Mordrate führt. Auch bei der Regression wird zwar häufig vom Einfluss von Variablen aufeinander gesprochen – doch auch hierbei handelt es sich um eine Zusammenhangsanalyse, die oftmals gleichzeitig erhobene Konstrukte betrachtet (z.B. aus Fragebogenstudien) und mit der Korrelation eng verwandt ist. So würden wir, wenn wir berechnen würden, ob COVID-Fälle durch 5G-Tower vorhergesagt werden, zu genau dem gleichen Ergebnis kommen, wie wenn wir die Vorhersagerichtung austauschen würden und die Anzahl der 5G-Tower durch COVID-Fallzahlen vorhersagen ließen. Auch die Regression kann per se keine Aussagen über Kausalzusammenhänge machen!
Um sich Rückschlüssen über die wahren Wirkungsbeziehungen aber zumindest anzunähern, können wir den Einfluss von Drittvariablen auf korrelative Zusammenhänge näher betrachten. (So wird z.B. die Korrelation zwischen 5G-Towern mit den COVID-Erkrankungen durch die Ballungsgebiete erklärt, in denen die Leute enger beieinander leben.) In dieser Sitzung beschäftigen wir uns daher mit der Partial- und der Semipartialkorrelation, d.h. Methoden, mit denen der Einfluss einer oder mehrerer Drittvariablen kontrolliert werden kann, um hierdurch Scheinkorrelationen, redundante oder maskierte Zusammenhänge aufzudecken.
Anmerkungen:
Vorbereitung des Datensatzes
Der Datensatz Depression
, den wir zur Illustration des heutigen Themas heranziehen werden, enthält fiktive Daten bezüglich Depressivitätswerten und einigen anderen, damit in Beziehung stehenden Variablen. Zusätzlich sind die kategorialen Variablen Intervention
und Geschlecht
enthalten, die wir in dieser Sitzung aber nicht weiter betrachten werden.
Variable | Kodierung |
---|---|
Lebenszufriedenheit | 1 - 10 |
Episodenanzahl | 1 - 10 |
Depressivitaet | 1 - 10 |
Neurotizismus | 1 - 10 |
Diesen Satz laden wir uns erst einmal herunter und prüfen die Struktur des Datensatzes:
load(url("https://pandar.netlify.app/daten/Depression.rda"))
str(Depression)
## 'data.frame': 90 obs. of 6 variables:
## $ Lebenszufriedenheit: num 7 5 8 6 6 8 7 5 8 6 ...
## $ Episodenanzahl : num 4 5 7 4 9 7 6 6 8 7 ...
## $ Depressivitaet : num 7 8 6 5 8 8 7 6 8 6 ...
## $ Neurotizismus : num 5 3 6 5 5 6 6 5 6 7 ...
## $ Intervention : Factor w/ 3 levels "Kontrollgruppe",..: 1 1 1 1 1 1 1 1 1 1 ...
## $ Geschlecht : Factor w/ 2 levels "0","1": 1 2 1 2 2 2 2 1 2 1 ...
Da wir die Variablen Intervention
und Geschlecht
für diese Sitzung nicht brauchen, entfernen wir diese mit dem uns bereits bekannten subset()
-Befehl:
dat <- subset(Depression, select = -c(Intervention, Geschlecht))
rm(Depression) # das nicht benötigte Objekt wird wieder gelöscht
Hier habe ich bei der Zuweisung in ein neues Objekt dat
als kürzeren Namen für unser Dataframe gewählt, um etwas Platz im Skript zu sparen und das Auge zu schonen. So können wir uns besser auf das Wesentliche konzentrieren.
Vorher müssen wir allerdings noch die Kodierung der Variable Neurotizismus
invertieren. Im ursprünglichen Datensatz wurde eine hohe Punktzahl als höhere emotionale Stablität (und damit geringerer Neurotizismus) kodiert, diese Kodierung wollen wir umkehren, damit eine höhere Punktzahl auch tatsächlich einen höheren Neurotizismus anzeigt:
# Neurotizismus invertieren
dat$Neurotizismus <- 11 - (dat$Neurotizismus)
Erinnerung: Korrelation
Wie Sie bereits aus Statistik 1 wissen, misst der Korrelationskoeffizient
Der Wert von
Empirisches Beispiel
Schauen wir uns zwei Variablen aus unserem Datensatz ab, bei denen die Vermutung nahe liegt, dass es einen (umgekehrten) Zusammenhang geben könnte: Lebenszufriedenheit und Depressivität. Sind Personen, die angeben, zufriedener mit ihrem Leben zu sein, weniger depressiv? Diese vermutete Korrelation lässt sich grafisch in einem Streudiagramm darstellen, was Sie mit dem Wissen aus den Lektionen zum ggplot2
-Paket wie folgt bewältigen können:
library(ggplot2) # notwendiges Paket für Grafiken laden
## Warning: Paket 'ggplot2' wurde unter R Version 4.3.2 erstellt
## Want to understand how all the pieces fit together? Read R
## for Data Science: https://r4ds.had.co.nz/
ggplot(dat, aes(x = Lebenszufriedenheit, y = Depressivitaet)) +
geom_point() +
theme_minimal() +
geom_smooth(method = "lm", se = FALSE)
## `geom_smooth()` using formula = 'y ~ x'

In dieser Abbildung mit einer dazugehörigen Regressionsgerade wird deutlich, dass hier deskriptiv ein negativer Zusammenhang vorliegt. Das heißt, eine höhere Lebenszufriedenheit geht mit einer geringeren Depressivität einher. Zu beachten ist hierbei natürlich, dass wir die x- und y-Achsen beliebig vertauschen könnten und aus diesem Befund daher keine Aussage über die Wirkungsrichtung ableitbar ist. Ob die Lebenszufriedenheit zu einer niedrigen Depressivität führt, oder eine niedrigere Depressivität nicht ursächlich für eine höhere Lebenszufriedenheit sein könnte, bleibt auf Basis dieser Daten komplett offen.
Mit cor.test()
können wir diesen Zusammenhang auch inferenzstatistisch überprüfen:
cor.test(dat$Depressivitaet, dat$Lebenszufriedenheit)
##
## Pearson's product-moment correlation
##
## data: dat$Depressivitaet and dat$Lebenszufriedenheit
## t = -5.0256, df = 88, p-value = 2.614e-06
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.6188046 -0.2938777
## sample estimates:
## cor
## -0.4722292
… und entnehmen dem Output entnehmen, dass die Korrelation
Die Beziehung zwischen der Stärke eines linearen Zusammenhangs und dem Anteil der gemeinsamen Varianz lässt sich durch Venn-Diagramme gut veranschaulichen:

Die Überlappung der beiden Kreise stellt die geteilte Varianz, d.h. die Korrelation dar, während die sich nicht überlappenden Flächen die voneinander unabhängig variierenden Anteile der Variablen ausmachen.
Drittvariablen
Der Zusammenhang zwischen zwei Variablen

In unserem Datensatz wollen wir uns als eine solche Drittvariable Neurotizismus
unter die Lupe nehmen. Während Lebenszufriedenheit und Depressivität in diesem Datensatz Zustände bezeichnen, die potenziell vorübergehend und zeitlich veränderlich sind (states). Im Gegensatz dazu ist Neurotizismus eine Persönlichkeitseigenschaft (trait), die eine stabile, langfristige Tendenz einer Person beschreibt, negative Emotionen wie Angst, Ärger und Traurigkeit zu erleben (vgl. Big-Five-Persönlichkeitsmodell nach Costa & McCrae, 1999).
Daher liegt es nahe, dass Neurotizismus im Zusammenhang mit Lebenszufriedenheit und Depressivität stehen könnte. Testen wir den vermuteten Zusammenhang in zwei einzelnen Korrelationstests:
# Korrelation der Drittvariablen mit den beiden ursprünglichen Variablen
cor.test(dat$Neurotizismus, dat$Lebenszufriedenheit)
##
## Pearson's product-moment correlation
##
## data: dat$Neurotizismus and dat$Lebenszufriedenheit
## t = -8.2276, df = 88, p-value = 1.587e-12
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.7623706 -0.5238208
## sample estimates:
## cor
## -0.659383
cor.test(dat$Neurotizismus, dat$Depressivitaet)
##
## Pearson's product-moment correlation
##
## data: dat$Neurotizismus and dat$Depressivitaet
## t = 12.289, df = 88, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.7035951 0.8603383
## sample estimates:
## cor
## 0.7948675
In der Tat zeigen die Ergebnisse, dass der Neurotizismus sowohl mit der Lebenszufriedenheit (-0.66) als auch mit Depressivität (0.79) statistisch bedeutsam korreliert.
Ist der Zusammenhang zwischen Lebenszufriedenheit und Depressivität also allein auf den Neurotizismus zurückzuführen? Oder besteht noch ein Zusammenhang zwischen den beiden, wenn man den Einfluss der Neurotizismus kontrolliert? Um dies zu überprüfen, müssen wir auf Methoden zur Kontrolle von Drittvariablen und zur Aufdeckung von Scheinkorrelationen, redundanten oder maskierten Zusammenhängen zurückgreifen: Die Partial- und Semipartialkorrelation.
Partialkorrelation
Die Partialkorrelation ist die bivariate Korrelation zweier Variablen

Denken wir nun an die Regression zurück, können wir die Partialkorrelation konzeptuell als etwas damit sehr verwandtes verstehen: Die Partialkorrelation
An unserem empirischen Beispiel werden wir das praktisch veranschaulichen.
Beispiel: Zusammenhang von Neurotizismus mit Lebenszufriedenheit und Depressivität
Nun heißt es näher zu betrachten, ob Neurotizismus eine mögliche konfundierende Variable in der Beziehung zwischen Lebenszufriedenheit und Depressivität darstellt. Diese Vermutung können wir mit einer Partialkorrelation überprüfen, für die wir die Residuen heranziehen. Diese bekommen wir aus einer einfachen Regression.
Erstmal stellen wir zwei einfache Regressionsmodelle mit Neurotizismus
als demjenigen Prädiktor auf, dessen Einfluss wir aus Depressivitaet
und Lebenszufriedenheit
isolieren möchten.
# Regression
reg_depr_neuro <- lm(Depressivitaet ~ Neurotizismus, data = dat)
reg_lz_neuro <- lm(Lebenszufriedenheit ~ Neurotizismus, data = dat)
Der Einfachheit halber legen wir die Residuen (also den von Neurotizismus
unabhängigen Varianzteil) in einem neuen Objekt ab:
# Residuen in Objekt ablegen (Anteil von Depr., der nicht durch Neuro. erklärt wird)
res_depr_neuro <- residuals(reg_depr_neuro)
# Residuen in Objekt ablegen (Anteil von LZ, der nicht durch Neuro. erklärt wird)
res_lz_neuro <- residuals(reg_lz_neuro)
Partialkorrelation als Korrelation der Residuen
Damit können wir die Partialkorrelation nun bestimmen:
cor(res_depr_neuro, res_lz_neuro)
## [1] 0.1137543
Mit einem Wert von Depressivitaet
und Lebenszufriedenheit
ohne Berücksichtigung der Drittvariablen – der ursprüngliche Zusammenhang verschwindet also nahezu.
Beachten Sie, dass wir hier NICHT cor.test()
verwenden, da die inferenzstatistische Absicherung die falschen Freiheitsgrade nutzen würde. Für interessierte Lesende findet sich in Appendix A eine genauere Erläuterung der Problematik.
Partialkorrelation mit Funktion
Natürlich müssen wir nicht jedes Mal die Residuen manuell berechnen, sondern können hierfür – nachdem Ihnen die Logik hinter der Partialkorrelation an diesem Beispiel deutlich werden konnte – auf eine Funktion zurückgreifen, die uns ppcor
installieren und das Paket laden müssen.
# Paket für Partial- und Semipartialkorrelation
install.packages("ppcor")
library(ppcor)
## Warning: Paket 'ppcor' wurde unter R Version 4.3.2 erstellt
## Lade nötiges Paket: MASS
##
## Attache Paket: 'MASS'
## Das folgende Objekt ist maskiert 'package:dplyr':
##
## select
Mit der Funktion pcor.test()
lässt sich die Partialkorrelation direkt ermitteln:
# Partialkorrelation mit Funktion
pcor.test(x = dat$Depressivitaet, # Das Outcome
y = dat$Lebenszufriedenheit, # Die Prädiktorvariable
z = dat$Neurotizismus) # wird aus X und Y auspartialisiert
## estimate p.value statistic n gp Method
## 1 0.1137543 0.2884924 1.067961 90 1 pearson
Vorsicht: Hier herrscht Verwechslungsgefahr bei der Notation der Argumente! Während wir typischerweise bei der Regression mit pcor.test
-Funktion in diesem Fall nichts damit zu tun. Mit x
ist hier die primäre Variable gemeint, aus deren gemeinsamer Varianz mit y
der Varianzanteil von z
herauspartialisiert werden soll. Das entspricht eher der Notation in den Venn-Diagrammen in dieser Lektion. Bei den Semipartialkorrelationen weiter unten wird diese Reihenfolge noch relevant!
In diesem Fall beträgt die Partialkorrelation
Unter Kontrolle des Neurotizismus verschwindet der ursprüngliche Zusammenhang zwischen Depressivität und Lebenszufriedenheit. Es handelt sich also um eine Scheinkorrelation zwischen den beiden Variablen, die durch den Neurotizismus erklärt wird. Anders ausgedrückt: Wie ausgeprägt die Depressivität einer Person ist, hängt nicht direkt mit ihrer Lebenszufriedenheit zusammen, sondern beide hängen damit zusammen, wie neurotisch die Person ist.
Übertragen auf die Regression: Wenn wir in ein einfaches Regressionsmodell mit Prädiktor Lebenszufriedenheit
zusätzlich den Prädiktor Neurotizismus
aufnehmen, wird das Regressionsgewicht Lebenszufriedenheit
nicht signifikant (
# Einfache lineare Regression von Depressivität auf Neurotizismus
mod1 <- lm(Depressivitaet ~ Lebenszufriedenheit, data = dat)
summary(mod1)$coefficients |> round(3)
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 9.122 0.709 12.862 0
## Lebenszufriedenheit -0.542 0.108 -5.026 0
# Erweiterte Regression, die Episodenanzahl einschließt
mod2 <- lm(Depressivitaet ~ Lebenszufriedenheit + Neurotizismus, data = dat)
summary(mod2)$coefficients |> round(3)
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.008 0.950 1.061 0.292
## Lebenszufriedenheit 0.105 0.099 1.068 0.288
## Neurotizismus 0.881 0.089 9.950 0.000
Auswirkung des Einschlusses von Drittvariablen
Wenn wir eine Drittvariable berücksichtigen und sich daraufhin die ursprüngliche Korrelation verändert, kann es zwei Szenarien geben: Entweder der Betrag der Korrelation (d.h. die Effektgröße ohne das Vorzeichen) wird kleiner, oder er wird größer. Hier finden Sie eine Interpretationshilfe, was dies bedeutet:
- Scheinkorrelation: Die Partialkorrelation wird (betraglich) kleiner als die ursprüngliche Korrelation (
)| r x y . z | < | r x y |
Wie in unserem Bespiel teilen alle drei Variablen miteinander Varianz. Partialisiert man nun eine Variable aus dem Zusammenhang der beiden anderen Variablen heraus, wird die geteilte Varianz weniger, womit die Korrelation betraglich sinkt (beachten Sie, dass es sich nur um die absolute Größe, d.h. den Betrag der Korrelation handelt, und nicht um die Vorzeichen!). Dieser Fall ist der am häufigsten eintretende, da in der Forschung oft Variablen auspartialisiert werden, weil es theoretische Annahme gibt, warum die Variablen Varianz teilen sollten, man aber eine isolierte Assoziation (Beziehung) zwischen
- Suppressoreffekt: Partialkorrelation ist (betraglich) größer als die ursprüngliche Korrelation (
)| r x y . z | > | r x y |
Es kann aber auch passieren, dass eine Drittvariable

Ein Beispiel: Sie untersuchen den Zusammenhang zwischen
Die habituelle Stimmung
Semipartialkorrelation
Es ist auch möglich, den Einfluss einer Drittvarible

Wenn wir unser Beispiel wieder aufgreifen, entspricht das konzeptuell der Korrelation zwischen Depressivität (ihrer gesamten Varianz) und dem Residualateil von Lebenszufriedenheit, der von Neurotizismus unabhängig ist:
# Semipartialkorrelation als Korrelation zwischen Depressivität (X) und
# Lebenszufriedenheit (Y), bereinigt um den Einfluss von Neurotizismus (Z) auf
# Lebenszufriedenheit (Y)
cor(dat$Depressivitaet, res_lz_neuro)
## [1] 0.06902416
Auch hier verwenden wir absichtlich die Funktion cor()
statt cor.test()
, da die inferenzstatistische Absicherung nicht korrekt durchgeführt werden würde aufgrund der Freiheitsgrade.
Mit der Funktion spcor.test()
aus dem zuvor installierten Paket ppcor
lässt sich die Semipartialkorrelation direkt ermitteln und die inferenzstatistische Absicherung gelingt.
# Semipartialkorrelation mit Funktion
spcor.test(x = dat$Depressivitaet, # Outcome
y = dat$Lebenszufriedenheit, # Prädiktor
z = dat$Neurotizismus) # wird aus Y, aber nicht X auspartialisiert
## estimate p.value statistic n gp Method
## 1 0.06902416 0.5203964 0.6453536 90 1 pearson
Der Koeffizient der Semipartialkorrelation
Wann brauche ich die Partial- und wann die Semipartialkorrelation?
Schließlich stellt sich natürlich die Frage, in welchem theoretischen Kontext man die Semipartialkorrelation gegenüber der Partialkorrelation überhaupt nutzen sollte.
Das kann sich in theoretischen Annahmen zur Wirkungsrichtung begründen. Bei der Partialkorrelation nehmen Sie an, dass die Drittvariable
Trotzdem gibt es noch einen Bereich, in dem sich die Semipartialkorrelation als nützlich erweist: im Kontext der Regression.
Rolle der Semipartialkorrelation in der Regression
Relevant wird das Konzept der Semipartialkorrelation im Rahmen der multiplen Regression bei der Berechnung des Determinationskoeffizienten
Erinnerung: Einfache Regression
In Statistik 1 haben wir bereits erfahren, dass in der einfachen linearen Regression der standardisierte
# Paket für standardisierte Beta-Koeffizienten
library(lm.beta)
## Warning: Paket 'lm.beta' wurde unter R Version 4.3.1 erstellt
# Einfache lineare Regression von Depressivität auf Neurotizismus
mod1 <- lm(Depressivitaet ~ Neurotizismus, data = dat)
lm.beta(mod1)$standardized.coefficients
## (Intercept) Neurotizismus
## NA 0.7948675
# Korrelation entspricht dem Beta-Koeffizienten
cor(dat$Depressivitaet, dat$Neurotizismus)
## [1] 0.7948675
Entsprechend ist die quadrierte Korrelation
Beispiel: Multiple Regression
Wir zeigen nun am Beispiel eines Modells mit zwei Prädiktoren – die in diesem Fall beide statistisch signifikanten Varianzanteil an der abhängigen Variable erklären – die Bedeutung der Semipartialkorrelation hierfür auf. Wir gehen hier wieder zur Notation über, in der
Bei multipler Regression ist der Determinationskoeffizient so definiert, dass für jeden weiteren Prädiktor das Quadrat einer Semipartialkorrelation aufaddiert wird, um den multiplen Determinationskoeffizienten zu bestimmen. Dabei wird der schon aufgenommene Prädiktor aus dem neuen herauspartialisiert und der Rest dann mit dem Kriterium korreliert (also die Semipartialkorrelation bestimmt):
In unserem Fall haben wir lediglich zwei Prädiktoren, wobei wir als Episodenanzahl
wählen. Für diese lassen wir uns direkt die standardisierten
# Ein Modell mit zwei Prädiktoren
mod2 <- lm(Depressivitaet ~ Neurotizismus + Episodenanzahl,
data = dat)
lm.beta(mod2) |> summary() # fügt std. Betas zum Output hinzu
##
## Call:
## lm(formula = Depressivitaet ~ Neurotizismus + Episodenanzahl,
## data = dat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.02852 -0.51230 -0.08009 0.65804 1.74145
##
## Coefficients:
## Estimate Standardized Std. Error t value Pr(>|t|)
## (Intercept) -1.17653 NA 0.46337 -2.539 0.0129
## Neurotizismus 0.82155 0.79784 0.05088 16.147 < 2e-16
## Episodenanzahl 0.51622 0.39471 0.06462 7.988 5.23e-12
##
## (Intercept) *
## Neurotizismus ***
## Episodenanzahl ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.8243 on 87 degrees of freedom
## Multiple R-squared: 0.7876, Adjusted R-squared: 0.7827
## F-statistic: 161.3 on 2 and 87 DF, p-value: < 2.2e-16
Wir legen den Determinationskoeffizienzen in ein Objekt ab:
R2_mod2 <- summary(mod2)$r.squared
…und berechnen auch die Korrelation zwischen corrs
(Korrelationen) ab.
r_yx1 <- cor(dat$Depressivitaet, dat$Neurotizismus)
r_yx2.x1 <- spcor.test(x = dat$Depressivitaet, # Outcome
y = dat$Episodenanzahl, # hier rauspartialisieren
z = dat$Neurotizismus)$estimate
corrs <- data.frame(r_yx1, r_yx2.x1)
corrs
## r_yx1 r_yx2.x1
## 1 0.7948675 0.3947031
Wie in der obigen Formel abstrakt dargestellt, entspricht auch in unserem konkreten Fall
sum(corrs^2) # Das entspricht dem Determinationskoeffizienten R^2
## [1] 0.7876049
R2_mod2
## [1] 0.7876049
Hierfür könnte man die Partialkorrelationen nicht verwenden – denn für jede Korrelation zwischen einem neuen hinzukommenden Prädiktor
Zusammenfassung
Wir haben in diesem Tutorial die Partial- und Semipartialkorrelation als Erweiterungen der Korrelation kennengelernt und ihre Beziehung zur multiplen Regression beleuchtet. Die Zusammenhänge zwischen den drei Größen soll der folgende Plot nochmal zusammenfassen.

Ausblick: Noch mehr Partialkorrelationen!
Das kann man natürlich nicht nur für zwei oder drei, sondern potentiell eine sehr große Anzahl von Variablen machen. Wenn wir in die Funktion pcor
mehrere Variablen übergeben, wird eine Matrix der Partialkorrelationen zurückgegeben. Hier ein Beispiel mit allen numerischen Variablen unseres Datensatzes. Die Funktion pcor()
gibt uns nicht nur die Tabelle der Partialkorrelationskoeffizienten aus, sondern auch die zugehörigen
# Tabelle der Partialkorrelationen
pcor_table <- pcor(dat)
pcor_table$estimate |> round(3)
## Lebenszufriedenheit Episodenanzahl
## Lebenszufriedenheit 1.000 -0.173
## Episodenanzahl -0.173 1.000
## Depressivitaet 0.199 0.662
## Neurotizismus -0.529 -0.565
## Depressivitaet Neurotizismus
## Lebenszufriedenheit 0.199 -0.529
## Episodenanzahl 0.662 -0.565
## Depressivitaet 1.000 0.825
## Neurotizismus 0.825 1.000
Man kann das ganze auch grafisch darstellen, um sich einen schnellen Überblick über ungerichtete Zusammenhänge zwischen einer größeren Anzahl von Variablen zu verschaffen. Dabei wird der Einfluss je aller anderen Variablen im Datensatz konstant gehalten (herauspartialisiert). Ein populäres Paket zur Visualisierung solcher ungerichteten Zusammenhänge ist das Paket qgraph
. Für die numerischen Variablen aus unserem Depressionsdatensatz könnte das Ganze so aussehen:
library(qgraph)
## Warning: Paket 'qgraph' wurde unter R Version 4.3.2 erstellt
# Partialkorrelationsnetzwerk grafisch darstellen:
Q <- qgraph(input = pcor_table$estimate,
layout = "spring", # grafische Anordnung der Variablen
edge.labels = TRUE, # Partialkorrelationen anzeigen
labels = substr(colnames(dat), 1, 5)) # die ersten 5 Zeichen der colnames

# Plot wird automatisch erzeugt, ohne dass wir Q ausführen müssen
Die Variablen werden automatisch durch die Einstellungen so angeordnet, dass Variablen, die stärker miteinander zusammenhängen, in der Abb. näher beieinander sind. Deskriptivstatistisch wird daran deutlich, dass z.B. der Zusammenhang zwischen Episodenlänge und Depressivität 0.66 beträgt, wenn der Einfluss von sowohl Neurotizismus als auch Lebenszufriedenheit herausgerechnet wird.
Ob die Zusammenhänge auch inferenzstatistisch abgesichert, d.h. auch in der Population von 0 verschieden sind, damit beschäftigen sich Forschende im Rahmen von Netzwerkanalysen. Diese Methode hat in der Persönlichkeits- und der klinischen Popularität in den letzten Jahren stark an Popularität gewonnen, um Zusammenhänge zwischen vielen simultan beobachteten, potenziell miteinander interagierenden Variablen zu untersuchen. Darüber können Sie in den pandaR-Lektionen zum KliPPs-Master noch mehr erfahren.
Appendix A
Inferenzstatistik der Partialkorrelation und weitere konzeptionelle Überlegungen
Bevor wir uns die Inferenzstatistik der Partialkorrelation genauer ansehen, wiederholen wir die Inferenzstatistik der Korrelation und schauen uns nochmal genau an, welche Eigenschaften eigentlich so ein Residuum hat.
Inferenzstatistik der Korrelation
Um eine Korrelation auf Signifikanz zu prüfen, untersuchen wir die Nullhypothese, dass die Korrelation
Liegt cor.test
Funktion werden genau mit dieser Formel, bzw. mit pt
(für den zugehörigen
# Infos aus cor.test
cortest <- cor.test(dat$Depressivitaet, dat$Lebenszufriedenheit)
# correlation r
r_cortest <- cortest$estimate
# t-Wert
t_cortest <- cortest$statistic
# p-Wert
p_cortest <- cortest$p.value
# zu Fuß:
r <- cor(dat$Depressivitaet, dat$Lebenszufriedenheit)
n <- nrow(dat) # Anzahl Personen
t <- r*sqrt(n-2)/sqrt(1-r^2)
p <- 2*pt(q = abs(t), df = n-2, lower.tail = F) # ungerichtet (deswegen *2)
# Vergleiche
# t
t_cortest
## t
## -5.025552
t
## [1] -5.025552
# p
p_cortest
## [1] 2.614365e-06
p
## [1] 2.614365e-06
Super, die Formeln stimmen also. Widmen wir uns nun den Regressionsresiduen.
Eigenschaften der Regressionsresiduen
Wir wollen
# Regression
reg_depr_neuro <- lm(Depressivitaet ~ Neurotizismus, data = dat) # x ~ z
reg_lz_neuro <- lm(Lebenszufriedenheit ~ Neurotizismus, data = dat) # y ~ z
# Residuen in Objekt ablegen (Residuen x)
res_depr_neuro <- residuals(reg_depr_neuro) # eps_x
# Residuen in Objekt ablegen (Residuen y)
res_lz_neuro <- residuals(reg_lz_neuro) # eps_y
Prüfen wir doch mal ein paar Eigenschaften des Residuums:
# Mittelwert 0
mean(res_depr_neuro) # mean(eps_x)
## [1] -9.714451e-17
round(mean(res_depr_neuro), 14) # mean(eps_x) gerundet auf 14 Nachkommastellen
## [1] 0
Das Residuum hat wie im Erwartungswert festgehalten den Mittelwert von 0. Weiterhin zeigt es auch keinen Zusammenhang zum Prädiktor in der Regression, in der es aufgetreten ist:
# Korrelation/Kovarianze mit dem Prädiktor z
cov(dat$Neurotizismus, res_depr_neuro) # 0
## [1] -7.573738e-17
cor(dat$Neurotizismus, res_depr_neuro) # 0
## [1] -4.109959e-17
round(cor(dat$Neurotizismus, res_depr_neuro), 16) # 0 gerundet auf 16 Nachkommastellen
## [1] 0
Nun interessiert uns aber, wie stark dieses Residuum noch mit der abhängigen Variable (
cor(dat$Depressivitaet, res_depr_neuro)
## [1] 0.6067831
Die Korrelation fällt relativ hoch aus:
Woher kommt die Formel der Partialkorrelation
Die Partialkorrelation hatten wir in diesem Beitrag über die Residuen einer Regression pro Variable definiert. Im vorherigen Abschnitt hatten wir die Eigenschaften dieser Residuen näher beleuchtet. Die Partialkorrelation war jetzt nichts weiter als die Korrelation der Residuen
# Partialkorrelation
partial_cor <- cor(res_lz_neuro, res_depr_neuro)
partial_cor
## [1] 0.1137543
Wir kennen aber auch die Formel für die Partialkorrelation:
Wenn wir die Korrelationen in Objekten abspeichern, so lässt sich diese Formel sehr leicht umsetzen:
# Partialkorrelation via Formel
r_xy <- cor(dat$Lebenszufriedenheit, dat$Depressivitaet)
r_xz <- cor(dat$Lebenszufriedenheit, dat$Neurotizismus)
r_yz <- cor(dat$Depressivitaet, dat$Neurotizismus)
r_xy.z <- (r_xy - r_xz*r_yz)/(sqrt(1-r_xz^2)*sqrt(1-r_yz^2))
r_xy.z
## [1] 0.1137543
Offensichtlich kommen beide Befehle zum gleichen Ergebnis. Interessierten sei gesagt, dass der verwendete Code für die Berechnung r_xy.z
das Gleiche ist wie:
r_xy.z <- (r_xy - r_xz*r_yz)/sqrt(1-r_xz^2)/sqrt(1-r_yz^2)
Hier besteht allerdings der Vorteil, dass man die Klammer im Nenner nicht vergessen kann.
Wir wollen uns nun der Frage widmen, wie diese Formel entsteht. Dazu müssen wir eine kleine Vorbereitung treffen: wir müssen den Depression
Datensatz standardisieren. Danach ist der Mittelwert und die Varianz jeder Variable 0 (Mittelwert) und 1 (Varianz). Dies geht sehr leicht mit dem scale
Befehl. Anschließend müssen wir alle Ananlysen nochmals wiederholen, aber keine Sorge, die Partialkorrelation bleibt identisch:
# Datensatz standardisieren:
dat_std <- data.frame(scale(dat))
# Regression
reg_depr_neuro <- lm(Depressivitaet ~ Neurotizismus, data = dat_std) # x ~ z
reg_lz_neuro <- lm(Lebenszufriedenheit ~ Neurotizismus, data = dat_std) # y ~ z
# Residuen in Objekt ablegen (Residuen x)
res_depr_neuro <- residuals(reg_depr_neuro) # eps_x
# Residuen in Objekt ablegen (Residuen y)
res_lz_neuro <- residuals(reg_lz_neuro) # eps_y
# Partialkorrelation
partial_cor <- cor(res_lz_neuro, res_depr_neuro)
partial_cor
## [1] 0.1137543
# Partialkorrelation via Formel
r_xy <- cor(dat_std$Lebenszufriedenheit, dat_std$Depressivitaet)
r_xz <- cor(dat_std$Lebenszufriedenheit, dat_std$Neurotizismus)
r_yz <- cor(dat_std$Depressivitaet, dat_std$Neurotizismus)
r_xy.z <- (r_xy - r_xz*r_yz)/(sqrt(1-r_xz^2)*sqrt(1-r_yz^2))
r_xy.z
## [1] 0.1137543
Sie sehen, alles bleibt identisch! Der Hauptunterschied ist nur, einige Korrelationen jetzt einfach Kovarianzen sind, da die Varianzen der Variablen = 1 sind:
cor(dat_std$Lebenszufriedenheit, dat_std$Depressivitaet)
## [1] -0.4722292
cov(dat_std$Lebenszufriedenheit, dat_std$Depressivitaet)
## [1] -0.4722292
Wir wissen, dass die Korrelation gerade die Kovarianz geteilt durch das Produkt der Standardabweichungen ist. Das bedeutet, dass in der Formel gelten muss:
Die Wurzel aus den letzten beiden Ausdrücken ergibt dann die Standardabweichung (von
# Cov[eps_x, eps_y]
cov(res_lz_neuro, res_depr_neuro)
## [1] 0.05189286
r_xy - r_xz*r_yz
## [1] 0.05189286
# -> identisch!
# Var[eps_x]
var(res_lz_neuro)
## [1] 0.5652141
1-r_xz^2
## [1] 0.5652141
# -> identisch!
# Var[eps_y]
var(res_depr_neuro)
## [1] 0.3681857
1-r_yz^2
## [1] 0.3681857
# -> identisch!
Daraus können wir nun ableiten, dass die Partialkorrelationsformel sich einfach aus der Korrelation der Residuen
Inferenzstatistik der Partialkorrelation
Nun kommen wir endlich zur Inferenzstatistik der Partialkorrelation. Dabei geht es hauptsächlich um die Frage, warum wir nicht einfach cor.test()
mit den Residuen nutzen können. Um eine Partialkorrelation auf Signifikanz zu prüfen, untersuchen wir die Nullhypothese, dass die Partialkorrelation
Liegt ppcor
Paket werden genau mit dieser Formel, bzw. mit pt
(für den zugehörigen
Wir sehen also, dass die cor.test
Funktion, wenn wir ihr die Residuen übergeben, gar nicht zum richtigen Ergebnis kommen kann,. Das liegt daran, dass die cor.test
Funktion davon ausgeht, dass ganz normale Variablen korreliert werden sollen. Residuen sind aber dahingehend besonders, dass sie bereits verrechnet wurden. Wir haben ja bereits Regressionen mit df = n-1-2
und nicht df=n-2
).
Weil es so viel Freude bereitet diese Inhalte tiefer zu verstehen, prüfen wir auch dies noch einmal:
# Infos aus pcor.test
pcortest <- pcor.test(dat$Depressivitaet, dat$Lebenszufriedenheit,
dat$Neurotizismus)
# correlation r
r_pcortest <- pcortest$estimate
# t-Wert
t_pcortest <- pcortest$statistic
# p-Wert
p_pcortest <- pcortest$p.value
# zu Fuß:
n <- nrow(dat) # Anzahl Personen
t <- r_xy.z*sqrt(n-1-2)/sqrt(1-r_xy.z^2)
p <- 2*pt(q = abs(t), df = n-1-2, lower.tail = F) # ungerichtet (deswegen *2)
# Vergleiche
# t
t_pcortest
## [1] 1.067961
t
## [1] 1.067961
# p
p_pcortest
## [1] 0.2884924
p
## [1] 0.2884924
Super, auch diesmal stimmen die Formeln. Im Laufe Ihres Studiums wird dieser Sachverhalt immer wieder aufgegriffen und auch noch näher erklärt. An dieser Stelle reicht es uns zu wissen, dass wir bei der Partialkorrelation andere