Trulli
from gourmet-blog.de

Überblick

Am Ende des Practicals wirst du wissen…

  1. Wie du einfache und multiple Regressionen in R implementierst.
  2. Wie du die Ergebnisse einer Regression interpretierst.

Aufgaben

A - Setup

  1. Öffne dein TheRBootcamp R project. Es sollte die Ordner 1_Data und 2_Code enthalten. Stelle sicher, dass du alle Datensätze, welche im Datensätze Tab aufgelistet sind, in deinem 1_Data Ordner hast.

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

## NAME
## DATUM
## Lineare Modelle Practical
  1. Speichere das neue Skript unter dem Namen lineare_modelle_practical.R im 2_Code Ordner.

  2. Lade das Paket tidyverse.

library(tidyverse)
  1. 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. Führe den Code unten aus um sicherzustellen, dass alle character Variablen als Faktoren vorliegen. Das hilft den statistischen Modellen kategoriale Variablen richtig zu interpretieren.

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

B - Einfache Regression

  1. In diesem ersten Abschnitt geht es darum den Effekt von Restzucker (Prädiktor) auf die wahrgenommene Qualität (Kriterium) des Weines zu prüfen. Verwende das template unten, um ein entsprechendes lineares Modell zu testen und in das Objekt wein_lm zu speichern.
# Einfache Regression
wein_lm <- lm(formula = XX ~ XX,
              data = XX)
wein_lm <- lm(formula = Qualität ~ Restzucker,
              data = wein)
  1. Printe das Fit-Objekt wein_lm. Welche Information wird dir angezeigt? Was sagen dir die Werte?
wein_lm

Call:
lm(formula = Qualität ~ Restzucker, data = wein)

Coefficients:
(Intercept)   Restzucker  
    5.85532     -0.00679  
  1. Der Default-Output von lm-Objekten ist nicht sehr informativ. Er zeigt nur den Code und die Schätzungen für die beiden Parameter: den Intercept und das Gewicht für Restzucker. Verwende summary() um einen informativeren Output zu erhalten.
# Zeige Ergebnisse
summary(XX)
summary(wein_lm)

Call:
lm(formula = Qualität ~ Restzucker, data = wein)

Residuals:
   Min     1Q Median     3Q    Max 
-2.851 -0.808  0.158  0.229  3.217 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  5.85532    0.01645  355.89   <2e-16 ***
Restzucker  -0.00679    0.00228   -2.98   0.0029 ** 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.873 on 6495 degrees of freedom
Multiple R-squared:  0.00137,   Adjusted R-squared:  0.00121 
F-statistic: 8.89 on 1 and 6495 DF,  p-value: 0.00287
  1. Summary zeigt einem ausführliche Ergebnisse in fünf Sektionen:
  • Call: Der Modellcode.
  • Residuals: Verschiedene Statistiken der Verteilung der Residuen.
  • Coefficients: Die Schätzer und Tests für die Modellparameter.
  • Signif. codes: Übliche Kodierung für Signifikanzlevels.
  • Schluss: Statistiken über dem Gesamtfit des Modells.

Welche Sektion zeigt an, welchen Einfluss Restzucker auf die Qualität des Weines hat?

  1. In der Tat, die Sektion Coefficients zeigt an, welche Werte für die Parameter geschätzt wurden und ob sie signifikant (gegeben H0) sind, und ist daher meist von zentralem Interesse. Welcher Paramterwert bzw. welches Gewicht (Estimate) wurde für Restzucker ermittelt? Was bedeutet dieser Wert?

  2. Exakt, ein Anstieg von 1 (g/ml) geht mit einer Veränderung von -.0067 in der Qualität einher. Also je süsser der Wein, desto niedriger die eingeschätzte Qualität. Aber ist dieser Effekt signifikant? Was zeigt der p-Wert für Restzucker an?

  3. Der Effekt von Restzucker ist signifikant von Null verschieden. Versuche nachzuverfolgen wie R zu diesem Wert kommt. Erinnere dich, der p-Wert folgt aus der Teststatistik, in diesem Fall t. Was ist der t-Wert für Restzucker?

  4. Der t-Wert ist -2.982. D.h., das beobachtete Ergebnis hat in der t-Verteilung mit 6,495 Freiheitsgraden (siehe Schlussabschnitt des Outputs) eine eben solche Ausprägung. Der p-Wert ist nun die Fläche unterhalb der Verteilung in die Richtung extremerer Werte. Die Funktion pt() gibt für bestimmte t-Werte die Fläche für niedrigere Werte raus. Probiert einmal was ihr erhaltet, wenn ihr den ermittelten t-Wert (q für Quantil) und die Freiheitsgrade (df) in die Funktion eingebt. Wie könnte der Zusammenhang zum p-Wert im Ergebnis aussehen?

# Berechne Fläche unter t-Wert
pt(q = XX, df = XX)
# Berechne Fläche unter t-Wert
pt(q = -2.982, df = 6495)
[1] 0.00144
  1. Genau der ermittelte Wert ist genau die Hälfte des p-Werts. Das hat damit zu tun, dass R-automatisch zweiseitig testet, d.h., Abweichungen in die positive und negative Richtung gleichermassen berücksichtigt. Hätte man eine klar gerichtete Hypothese, könnte man also den p-Wert im Ergebnis halbieren.

  2. Es gibt noch weitere Werte, die du miteinander verknüpfen kannst. Der t-Wert ergibt sich nämlich unmittelbar aus dem Gewicht (Estimate) und dem Standardfehler (Std.Error). Kannst du den Zusammenhang erkennen?

  3. Der t-Wert ist nichts anderes als das Gewicht geteilt durch den Standardfehler. Betrachte nun die andere, erste Zeile des Coefficients Outputs. Wie interpretierst du ihn? Die Zeile gibt die Ergebnisse für den Intercept an, d.h., den Wert des Kriterium, wenn der Prädiktor Null ist. Nicht leicht zu interpretieren, oder? Es gibt keinen Wein mit Restzucker gleich Null. Später mehr dazu.

  4. Lenke nun deine Aufmerksamkeit auf den Schlussteil. Hier siehst du eine weitere Zusammenfassung über das Modell. Von besonderem Interesse ist das Multiple R-squared, da es angibt wie viel Varianz das Modell insgesamt am Kriterium erklären kann. Was meinst du, ist der Wert gross oder klein?

  5. Der Wert ist in der Tat relativ klein. Nur .1% der Varianz kann aufgeklärt werden, was wirklich nicht viel ist. Zeit weitere Prädiktoren hinzuzunehmen.

C - Multiple Regression

  1. Teste nun den Effekt mehrerer Prädiktoren auf die wahrgenommene Qualität (Kriterium) des Weines. Inkludiere neben Restzucker pH_Wert, Alkohol, und Sulphate. Verwende das template unten, um ein entsprechendes lineares Modell zu testen und in das Objekt wein_lm zu speichern.
# Multiple Regression
wein_lm <- lm(formula = XX ~ XX + XX + XX + XX,
              data = XX)
wein_lm <- lm(formula = Qualität ~ Restzucker + pH_Wert + Alkohol + Sulphate,
              data = wein)
  1. Printe wiederum das Modell. Wie interpretierst du die Grössen der Regressionsgewichte Kannst du aus ihnen ablesen, welcher Prädiktor die meiste Bedeutung hat?
wein_lm

Call:
lm(formula = Qualität ~ Restzucker + pH_Wert + Alkohol + Sulphate, 
    data = wein)

Coefficients:
(Intercept)   Restzucker      pH_Wert      Alkohol     Sulphate  
     1.8422       0.0280      -0.0767       0.3669       0.4171  
  1. Rohe Regressionsgewichte lassen nicht darauf schliessen, welche Bedeutung ein Prädiktor hat, da die Skalierung der Prädiktoren unmittelbar mit einfliesst. Wäre z.B. Alkohol in Promille angegeben, dann wäre das Gewicht 10x grösser. Dies kannst du leicht mit dem Template unten testen. Füge die vier Prädiktoren ein sodass Alkohol in innerhalb der I() steht und führe das Modell aus. Speichere es aber nicht.
# Multiple Regression
lm(formula = XX ~ XX + XX + I(XX / 10) + XX,
   data = XX)
lm(formula = Qualität ~ Restzucker + pH_Wert + I(Alkohol/10) + Sulphate,
   data = wein)

Call:
lm(formula = Qualität ~ Restzucker + pH_Wert + I(Alkohol/10) + 
    Sulphate, data = wein)

Coefficients:
  (Intercept)     Restzucker        pH_Wert  I(Alkohol/10)       Sulphate  
       1.8422         0.0280        -0.0767         3.6687         0.4171  
  1. Also, hohe Regressionsgewichte müssen immer relativ zur Skalierung interpretiert werden. Nichtsdestotrotz ist bereits etwas bemerkenswertes an den Gewichten ablesbar. Schau dir das Regressionsgewicht für Restzucker an. War das nicht mal negativ? Führt mehr Zucker plötzlich zu höherer wahrgenommener Qualität?

  2. So einfach ist das leider nicht. In der multiplen Regression hängt alles mit allem zusammen. D.h., das Gewicht von Restzucker hängt von den anderen Gewichten mit ab und ist konditional zu interpretieren. Wiederum nicht leicht!

  3. Schaue dir nun die summary() an. Welche Effekte gehen in welche Richtung, welche Effekte sind signifikant?

summary(wein_lm)

Call:
lm(formula = Qualität ~ Restzucker + pH_Wert + Alkohol + Sulphate, 
    data = wein)

Residuals:
   Min     1Q Median     3Q    Max 
-3.494 -0.501 -0.035  0.492  3.099 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  1.84215    0.22155    8.32  < 2e-16 ***
Restzucker   0.02800    0.00225   12.43  < 2e-16 ***
pH_Wert     -0.07671    0.06256   -1.23     0.22    
Alkohol      0.36687    0.00864   42.49  < 2e-16 ***
Sulphate     0.41713    0.06646    6.28  3.7e-10 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.772 on 6492 degrees of freedom
Multiple R-squared:  0.219, Adjusted R-squared:  0.219 
F-statistic:  456 on 4 and 6492 DF,  p-value: <2e-16
  1. Alkohol, Sulphate, und Restzucker hängen alle drei signifikant positiv mit der Qualität zusammen, der pH-Wert dagegen negativ, jedoch nicht signifikant. Auf Basis der Summary, kannst du jetzt ausmachen, welcher Prädiktor die grösste Bedeutung hat? Welcher Wert würdest du zu Rate ziehen?

  2. Am besten zieht man in der Tat den t-Wert heran. Je grösser dieser ist, desto grösser ist das Signal für den assoziierten Prädiktor. Richte nun deine Aufmerksamkeit auf den Schlussteil. Wie gross ist R-squared?

  3. Das R-squared hat sich deutlich erhöht. Nun werden knapp 22% der Varianz des Kriteriums erklärt. Auch dies ist etwas, dass wir unmittelbar mit den Funktionen nachvollziehen können. Speichert mit dem template zunächst das Kriterium, die Residuen und die gefitteten Werte. Für letztere beiden gibt es, wie ihr seht, nützliche Funktionen.

# Residuen & gefittete Werte
kriterium <- wein$XX
residuen <- residuals(XX)
gefitted <- fitted(XX)
# Residuen & gefittete Werte
kriterium <- wein$Qualität
residuen <- residuals(wein_lm)
gefitted <- fitted(wein_lm)
  1. Jetzt teilt ihr einmal die Varianz der gefitteten Werte (var(gefitted)) durch die Varianz des Kriteriums (var(kriterium)).

  2. Et voilà, ihr habt R-squared gefunden. Ihr könnt jetzt auch einmal die Varianz der Residuen und die der gefitteten Werte addieren und mit der Varianz des Kriteriums vergleichen.

D - Multiple Regression: Alle Prädiktoren

  1. Teste nun den Effekt aller Prädiktoren (bis auf Farbe) auf die wahrgenommene Qualität (Kriterium) des Weines. Um alle Prädiktoren zu inkludieren kann man einfach . verwenden, und um einen auszuschliessen -- Verwende das template unten, um ein entsprechendes lineares Modell zu testen und in das Objekt wein_lm zu speichern.
# Multiple Regression aller Prädiktoren ausser Farbe
wein_lm <- lm(formula = XX ~ . - XX,
              data = XX)
# Multiple Regression aller Prädiktoren ausser Farbe
wein_lm <- lm(formula = Qualität ~ . - Farbe,
              data = wein)
  1. Evaluiert die Summary des Modells.
  • Sind die Vorzeichen der Gewichte gleichgeblieben?
  • Welche Prädiktoren sind signifikant?
  • Welche Prädiktoren haben den grössten Einfluss?
  • Um wie viel ist R-squared gestiegen?

E - Multiple Regression: Alle Prädiktoren

  1. Gehe zurück zum ersten Modell mit Restzucker als einzigen Prädiktor. Inkludiere nacheinander jeweils einen anderen Prädiktor und studiere, wie sich das Gewicht und die Signifikanz von Restzucker verändert.
# Multiple Regression aller Prädiktoren ausser Farbe
wein_lm <- lm(formula = XX ~ Restzucker + XX,
              data = XX)
# Multiple Regression aller Prädiktoren ausser Farbe
lm(formula = Qualität ~ Restzucker + Alkohol,
   data = wein)

Call:
lm(formula = Qualität ~ Restzucker + Alkohol, data = wein)

Coefficients:
(Intercept)   Restzucker      Alkohol  
     1.8754       0.0259       0.3624  
  1. Der Schlüssel zur letzten Aufgabe steckt in cor(wein %>% select(-Farbe))['Restzucker',] und in cor(wein %>% select(-Farbe))['Qualität',], d.h., den bivariaten Korrelationen zwischen Restzucker und den jeweils anderen Prädiktoren und Qualität und den jeweils anderen Prädiktoren. Beides fliesst in die Berechnung des Gewichts für Restzucker mit ein. Siehe Folien.

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