Lineare Modelle I

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.

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

3. Speichere das neue Skript unter dem Namen lineare_modelle_practical.R im 2_Code Ordner.

4. Lade das Paket tidyverse.

library(tidyverse)

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

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

6. Printe den Datensatz.

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

8. 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 rechnen 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)

2. Printe das Fit-Objekt wein_lm. Welche Informationen werden dir angezeigt? Was sagen dir die Werte?

wein_lm

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

Coefficients:
(Intercept)   Restzucker  
    5.85532     -0.00679  

3. Der Default-Output von lm-Objekten ist nicht sehr informativ. Er zeigt nur den Code und die Schätzungen für die zwei 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

4. 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?

5. Die Sektion Coefficients zeigt an, welche Werte für die Parameter geschätzt wurden und ob sie signifikant (gegeben H0) sind. Welcher Paramterwert bzw. welches Gewicht (Estimate) wurde für Restzucker ermittelt? Was bedeutet dieser Wert?

6. 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?

7. Der Effekt von Restzucker ist signifikant von Null verschieden. Beginne 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?

8. 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 der relative Anteil möglicher, extremerer Werte in der Verteilung. Die Funktion pt() gibt für bestimmte t-Werte die Anteile der niedrigeren 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

9. 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.

10. 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?

11. Der t-Wert ist nichts anderes als das Gewicht geteilt durch den Standardfehler. Betrachte nun die erste Zeile des Coefficients Outputs. Wie interpretierst du sie? 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 nämlich keinen Wein mit Restzucker gleich Null. Später mehr dazu.

12. Lenke nun deine Aufmerksamkeit auf den Schlussteil. Hier siehst du eine Zusammenfassung über die Gesamtqualität des Modells. Von besonderem Interesse ist das Multiple R-squared. Es gibt an wie viel Varianz das Modell insgesamt am Kriterium erklären kann. Was meinst du, ist der Wert gross oder klein?

13. Der Wert ist eher 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 auch pH_Wert, Alkohol, und Sulphate. Verwende das Template unten, um ein entsprechendes lineares Modell zu rechnens und al 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)

2. 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  

3. 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 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  

4. Rgressionsgewichte müssen immer relativ zur Skalierung interpretiert werden. Nichtsdestotrotz ist bereits etwas interessantes 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?

5. 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, was bei mehreren Prädiktoren schwierig werden kann!

6. 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

7. 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? Welchen Wert würdest du zu Rate ziehen?

8. Am besten zieht man 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?

9. 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 nachvollziehen können. Speichert mit dem Template unten 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)

10. Jetzt teilt ihr einmal die Varianz der gefitteten Werte (var(gefitted)) durch die Varianz des Kriteriums (var(kriterium)).

11. 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 rechnen und als 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)

2. 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: Gewichte verstehen

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  

2. 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 die 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

• Discovering Statistics with R von Andy Field ist sehr gut
• YaRrr! The Pirate’s Guide to R hat hilfreiche und unterhaltsame Kapitel zu Statistik mit R.