Trulli
from delinat.com

Überblick

Am Ende des Practicals wirst du wissen…

  1. Wie du kategoriale Variablen als Prädiktoren in einer Regression verwendest und interpretierst.
  2. Wie du Interaktionen integrierst und warum es wichtig sein kann die Variablen zu standardisieren.

Aufgaben

A - Setup

  1. Öffne dein TheRBootcamp R project.

  2. Öffne ein neues R Skript. Schreibe deinen Namen, das Datum und “Lineare Modelle II Practical” als Kommentare an den Anfang des Skripts.

## NAME
## DATUM
## Lineare Modelle II Practical

  1. Speichere das neue Skript unter dem Namen lineare_modelle_II_practical.R im 2_Code Ordner.

  2. Lade das Paket tidyverse.

  3. Verwende die read_csv() Funktion um wein.csv einzulesen.

# Lese Daten ein
wein <- read_csv(file = "1_Data/wein.csv")

  1. Printe den Datensatz.

  2. Verwende names(XX), summary(XX), und View(XX) um einen weiteren Überblick über die Daten zu bekommen.

  3. Wiederum, führe den Code unten aus um sicherzustellen, dass alle character Variablen als Faktoren vorliegen, was den statistischen Modellen hilft kategoriale Variablen richtig zu interpretieren.

# Konvertiere alle character zu factor
wein <- wein %>% mutate_if(is.character, factor)

B - Gruppenvergleiche: t-Test versus Regression

  1. In diesem Abschnitt geht es darum den Effekt der Farbe (Prädiktor) des Weins, rot oder weiss, auf die Qualität des Weines zu prüfen. Verwende hierzu zunächst den Code unten um zwei Vektoren zu erstellen, die jeweils die Qualitätswerte für weisse und rote Weine beinhalten.
# Qualitätsvektoren nach Farbe
weiss <- wein %>% filter(Farbe == 'weiss') %>% pull(Qualität)
rot <- wein %>% filter(Farbe == 'rot') %>% pull(Qualität)
  1. Verwende nun das t.test() template unten um die beiden Vektoren mit einem t-Test miteinander zu vergleichen. Du brauchst das Ergebnis nicht zu speichern.
# t-test
t.test(x = XX, y = XX)

  1. Was verrät dir der Output des t-Tests über den Unterschied weisser und roter Weine in der wahrgenommenen Qualität? Die Antwort findest du in der zweiten Zeile (beginnend mit t=..) zusammen mit der letzten Zeile.

  2. Weisse Weine wurden 0.241886 (Differenz der beiden Mittelwerte) Punkte als qualitativ wertiger wahrgenommen und dieser Unterschied ist signifikant. Versuche nun das selbe Ergebnis mit einer Regression zu erzielen. Prädiziere hierzu Qualität durch Farbe.

# Regression
wein_lm <- lm(formula = XX ~ XX, 
              data = XX)

  1. Printe zunächst einmal das Objekt und betrachte die Regressionsgewichte. Kommen dir eine oder sogar beide der Zahlen bekannt vor?

  2. Genau, das Regressionsgewicht für Farbe spiegelt exakt den Unterschied der beiden Mittelwerte wieder. Was repräsentiert dann der Intercept?: Den Wert derjenigen Kategorie, der R den Wert 0 zugewiesen hat. Dies ist bei ursprünglich character Variablen per Default diejenige Kategorie die alphabetisch früher kommt, d.h., 'rot' < 'weiss'.

  3. Lasse dir jetzt die summary() für das Modell ausgeben und vergleiche Freiheitsgrade, t- und p-Wert mit denen aus dem t-Test.

  4. Die Werte im t-Test stimmen trotz der gleichen Schätzwerte für die Gewichte nicht ganz überein. Dies liegt daran, dass die t.test() Funktion per Default eine Variante rechnet die erlaubt, dass die Varianzen in den beiden Gruppen (rot und weiss) unterschiedlich sind. Die Regression dagegen geht stets davon aus, dass diese gleich sind, was wir im Detail in der Session zu robuste Statistiken besprechen werden. Rechne den t-Test noch einmal, diesmal mit dem zusätzlichen Argument var.equal = TRUE.

# t-test
t.test(x = XX, y = XX, var.equal = )
  1. Ordnung ist eingekehrt. Alle Werte im t-test und der Regression sollten jetzt identisch sein, was daran liegt, dass die beiden Modelle nun exakt identisch sind.

C - Gruppenvergleiche: Kodierung

  1. Per Default verwendet R die Dummy-Kodierung, in der eine Ausprägung den Wert 0 erhält und eine andere den Wert 1. Eine alternative Kodierung existiert in der Effektkodierung, die stattdessen die Werte -1 und 1 vergibt. Um die Konsequenzen dieser beiden Kodierungen direkt zu vergleichen, erstelle zunächst für jede Kodierung eine neue Variable im Datensatz unter Verwendung des Codes unten.
# Kodierungen der Farbe
wein <- wein %>% mutate(Farbe_dummy = ifelse(Farbe == 'rot', 0, 1),
                        Farbe_effekt = ifelse(Farbe == 'rot', -1, 1))
  1. Rechne nun zwei Regressionen, jeweils eine mit einer der beiden Farbkodierungen als Prädiktor und speichere sie ab als wein_dummy und wein_effekt.
# Regression dummy
wein_dummy <- lm(formula = XX ~ XX, 
                 data = XX)

# Regression effekt
wein_effekt <- lm(formula = XX ~ XX, 
                  data = XX)

  1. Printe nun die beiden Objekte und vergleiche die Gewichte. Diejenigen der Dummy-Kodierung sollten dir bekannt vorkommen. Wie verhalten sich die Gewichte der Effektkodierung dagegen? Erkennst du den Zusammenhang?

  2. Der Effekt der Farbe unter der Effektkodierung entspricht genau der Hälfte des Effekts unter der Dummykodierung und um dies auszugleichen hat sich der Intercept auch um genau den selben Wert verändert, nur in die andere Richtung. Verifiziere die Gewichte mit den einfachen Rechnungen im Code unten.

# Dummy-Kodierung
mean(rot) # intercept
mean(weiss) - mean(rot) # gewicht farbe

# EffekKodierung
(mean(rot) + mean(weiss))/2 # intercept
mean(weiss) - (mean(rot) + mean(weiss))/2 # gewicht farbe

  1. Vergleiche nun die summary() der beiden Modelle. Was unterscheidet sich, und was nicht?

  2. Die Kodierung hat ausgenommen die Skalierung des Regressionsgewichts und des zugehörigen Standardfehlers keinen Einfluss. t-Wert und Signifikanz sind identisch. Einzig der Intercept ändert sich wesentlich, da er unter der Effektkodierung tatsächlich weiter von Null entfernt ist.

D - Interaktionen

  1. Der erste Abschnitt hat gezeigt, dass weisser Wein gegenüber Rotem bevorzugt wurde. Könnte dies daran liegen, dass nicht für andere Variablen kontrolliert wurde? Oder könnte es sein, dass dieser Unterschied nur unter bestimmten Bedingungen gilt, d.h., dass es Moderatoren für den Effekt der Farbe gibt? Rechne eine Regression die zusätzlich Alkohol mit aufnimmt und verknüpfe die Prädiktoren mit * sodass auch die Interaktion berücksichtigt wird.
# Regression mit Interaktion
wein_lm <- lm(formula = XX ~ XX * XX, 
              data = XX)

  1. Printe das Objekt und betrachte die Gewichte. Wie interpretierst du die einzelnen Werte?

  2. Die Gewichte lassen sich folgendermassen interpretieren: Unter Berücksichtigung von Alkohol und der Interaktion, werden weisse Weine um .7 besser bewertet. Unter Berücksichtigung der Farbe (d.h. für Farbe == ‘rot’) und der Interaktion, führt ein Anstieg von einem Volumenprozent zu einer Verbesserung der wahrgenommenen Qualität von .36. Unter Berücksichtigung von Farbe und Alkohol, führt die Interaktion, also das Produkt der beiden Prädiktoren zu einer Veränderung von -.05. Dies ist so zu interpretieren, dass der Effekt von Alkohol unter Weissweinen um genau diese Grösse kleiner ist. Dies bedeutet im Umkehrschluss dass der Effekt von Farbe unter hoch-alkoholischen Weinen kleiner ausfällt. Verwendet den Code unten um diese Ergebnisse zu visualisieren. Gelb bedeutet Weisswein, lila Rotwein.

# Visualisierung
ggplot(data = wein, 
       aes(x = Alkohol, y = Qualität, col = Farbe, alpha = .01)) + 
  scale_color_manual(values = viridis::cividis(2)) + 
  geom_jitter(width=2,height=1.5,size=.1) + theme_minimal() + theme(legend.position = 'none') +
  geom_smooth(data = wein %>% filter(Farbe == 'weiss'), method = 'lm') + 
  geom_smooth(data = wein %>% filter(Farbe == 'rot'), method = 'lm')

  1. Lasse dir nun die summary() anzeigen und inspiziere die t-Werte und Signifikanzen. Welche Prädiktoren sind signifikant?

  2. Alle drei Prädiktoren sind signifikant. t-Wert und Signifikanz sind am extremsten für Alkohol, jedoch ist das Gewicht für Alkohol nur halb so gross. Wie ist dies zu erklären?

  3. Genau die Gewichte hängen auch von der Skalierung der Prädiktoren ab. Rechne nun die Regression noch einmal, aber diesmal mit skalierten Prädiktoren. Siehe Code unten.

# Skalierungsfunktion
skaliere = function(x) (x - mean(x))/sd(x)

# Regression mit skalierten Prädiktoren
wein_lm <- lm(formula = XX ~ XX * XX, 
              data = XX %>% mutate_if(is.numeric, skaliere))

  1. Printe zunächst einmal das Objekt und betrachte die Gewichte. Alle Werte haben sich substantiell verändert. Insbesondere hat nun Alkohol das zweifach höhere Gewicht. Wie sind diese neuen Gewichte zu interpretieren, wenn du bedenkst, dass nach der Skalierung alle Variablen eine Standardabweichung = 1 haben?

  2. Nach Skalierung sind die Gewichte als Veränderungen in Standardabweichungen zu interpretieren. Z.B., eine Veränderung von einer Standardabweichungen in Alkohol führt demnach zu einer Veränderung von .4928 Standardabweichungen in der Qualität. Am meisten hat sich der Intercept verändert. Wie interpretiert ihr diesen unter der Skalierung?

  3. Der Intercept ist der geschätzte Mittelwert des Kriteriums für den Fall in dem alle Prädiktoren Null sind. Ohne Skalierung ist dieser Fall schwer zu interpretieren, da Null nicht vorkommen muss in den Prädiktoren. Unter der Skalierung ändert sich dies aber, denn nun haben alle Prädiktoren, mit Ausnahme des Faktors Farbe einen Mittelwert von Null. Das heisst der Intercept spiegelt die Qualität für ein mittleres Level von Alkohol, rote Weine und keine Interaktion in Standardabweichungen wieder. Dies kannst du wiederum in der Visualisierung nachvollziehen. Benutze den Code unten und suche den Punkt auf der lila Linie für Alkohol = 0.

# Visualisierung
ggplot(data = wein %>% mutate_if(is.numeric, skaliere), 
       aes(x = Alkohol, y = Qualität, col = Farbe, alpha = .01)) + 
  scale_color_manual(values = viridis::cividis(2))  +
  geom_jitter(width=2,height=1.5,size=.1) + theme_minimal() + theme(legend.position = 'none') +
  geom_smooth(data = wein %>% mutate_if(is.numeric, skaliere) %>% filter(Farbe == 'weiss'), method = 'lm') + 
  geom_smooth(data = wein %>% mutate_if(is.numeric, skaliere) %>% filter(Farbe == 'rot'), method = 'lm')

  1. Schaue dir nun die summary() an. Was hat sich verändert bzgl. t-Werten, Signifikanzen, und Standardfehlern, und R-squared?

  2. Ein wesentlicher Einfluss besteht auf die Variable Farbe, welche als einzige nicht mit standardisiert wurde. Vor der Standardisierung von Alkohol war Farbe beinahe perfekt mit der Interaktion korreliert, was zu einem extrem aufgeblasenen Standardfehler geführt hat. Die Standardisierung hat diese Korrelation jedoch stark reduziert, so dass der Standardfehler auf das Niveau der anderen Prädiktoren gefallen ist. In der Konklusion lohnt es sich beinahe immer die Prädiktoren (oder auch das Kriterium) zu standardisieren, wenn ein Vergleich der Gewichte und Signifikanzen von Interesse ist.

Beispiele

# Regression mit R

library(tidyverse)

# Model:
# Sagt der Hubraum (displ) die pro gallone 
# fahrbaren Meilen voraus?
hwy_mod <- lm(formula = hwy ~ displ,
               data = mpg)

# Ergebnisse 
summary(hwy_mod)
coef(hwy_mod)

# Gefittete Werte
hwy_fit <- fitted(hwy_mod)
hwy_fit

# Residuums 
hwy_resid <- residuals(hwy_mod)
hwy_resid

Datensätze

Datei Zeile Spalte
wein.csv 6497 13

wein.csv

Der wein.csv Datensatz enthält aus den Jahren 2004 bis 2007 des Comissão De Viticultura Da Região Dos Vinhos Verdes, der Offiziellen Zertifizierungsagentur des Vinho Verde in Portugal.

Name Beschreibung
Qualität Qualitätsurteil über den Wein von 1-9
Farbe Roter oder weisser Wein
Gelöste_Säure Konzentration der im Wein gelösten Säuren
Freie_Säure Konzentration der verflüchtigbaren Säuren
Citronensäure Citronensäurekonzentration im Wein
Restzucker Zuckerkonzentration im Wein
Chloride Chloridkonzentration im Wein
Freie_Schwefeldioxide Konzentration der verflüchtigbaren Schwefeldioxide
Gesamt_Schwefeldioxide Konzentration der Schwefeldioxide insgesamt
Dichte Dichte des Weins
pH_Wert pH-Wert des Weins. Je kleiner, desto saurer.
Sulphate Sulphatkontration im Wein
Alkohol Alkoholkonzentration im Wein in %

Funktionen

Pakete

Package Installation
tidyverse install.packages("tidyverse")

Funktionen

Function Package Description
lm stats Fitte ein lineares Modell
fitted stats Extrahiere vorhergesagte Werte
residuals stats Extrahiere Residuen

Resourcen

Books