Statistics with R The R Bootcamp |
Am Ende des Practicals wirst du wissen…
Öffne dein TheRBootcamp
R project.
Öffne ein neues R Skript. Schreibe deinen Namen, das Datum und “Lineare Modelle III Practical” als Kommentare an den Anfang des Skripts.
## NAME
## DATUM
## Lineare Modelle III Practical
Speichere das neue Skript unter dem Namen lineare_modelle_III_practical.R
im 2_Code
Ordner.
Lade die Pakete tidyverse
und MASS
.
Verwende die read_csv()
Funktion um wein.csv
und bigmart.csv
einzulesen.
# Lese Daten ein
wein <- read_csv(file = "1_Data/wein.csv")
bigmart <- read_csv(file = "1_Data/bigmart.csv")
Printe die Datensätze.
Verwende names(XX)
, summary(XX)
, und View(XX)
um einen weiteren Überblick über die Daten zu bekommen.
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)
bigmart <- bigmart %>% mutate_if(is.character, factor)
Farbe
vorhersagen. In anderen Worten, in welchen Gesichtspunkten unterscheiden sich weisse und rote Weine voneinander. Verwende das Template um in einer Regression Farbe
durch Alkohol
, Dichte
, und Sulphate
vorherzusagen. Achtung: Es tritt ein Fehler auf. Siehe nächste Aufgabe.# Regression
wein_lm <- lm(formula = XX ~ XX + XX + XX,
data = wein)
'binomial'
Link und Modell zu fitten, aka eine logistsiche Regression. Nenne das entstehende Objekt wein_glm
.# Logistische Regression
wein_glm <- glm(formula = XX ~ XX + XX + XX,
data = XX,
family = 'XX')
Printe wein_glm
. Was verrät dir der Output?
Zu allererst ist dir vielleicht aufgefallen, dass einige der Gewichte sehr gross sind. Vor allem der Intercept ist erstaunlich hoch, wenn man bedenkt dass Farbe eigentlich nur Werte zwischen 0 und 1 annehmen kann. Dies hat natürlich mit dem Link-Funktion zu tun. Das bedeutet, im generalisierten linearen Modell muss nun auch noch die Link-Funktion mitberücksichtigt werden.
Weiterhin ist die sicher aufgefallen, dass am Schluss des Outputs Ergebnisse angezeigt werden, die bei einem einfachen linearen Modell fehlen. Diese Werte haben etwas damit zu, dass nun mit Maximum Likelihood gefittet wird. Später mehr dazu.
Lasse dir nun die summary()
anzeigen. Der Output ist sehr ähnlich aufgebaut wie der von lm()
. Lenke deine Aufmerksamkeit zunächst auf den Teil zu Coefficients
. Welche Prädiktoren sind signifikant?
Alle Prädiktoren tragen signifikant zur Vorhersage der Farbe
bei. Beachte, dass hier z-Werte anstatt t-Werte angegeben sind, was darauf zurückzuführen ist, dass hier ein anderes Datenmodell verwendet wird. Schaue dir nun den Schlussteil des Outputs an. Analog zur Standard Regression stehen hier Informationen über den Gesamtfit des Modells. Zunächst einmal wird die Null deviance
angegeben. Dies ist die Performanz des sogenannten Null-Modells, welches nur einen Intercept verwendet. Dieses Modell sagt für jeden Wein schlicht die relative Häufigkeit der häufigeren Kategorie voraus. D.h., im Datensatz ist der Anteil weisser 0.754. Entsprechend sagt ein solches Modell für jeden Wein voraus, dass eine Wahrscheinlichkeit von 0.754, dass der Wein weiss ist. Der folgende Code zeigt, wie du damit auf die Null deviance
kommst.
# Wahrscheinlichkeit weiss für Null-Modell
wsk_weiss <- mean(wein$Farbe == 'weiss')
# Log likelihood Null-Modell
log_likelihood <- sum(log(wsk_weiss) * sum(wein$Farbe == 'weiss')) +
sum(log(1-wsk_weiss) * sum(wein$Farbe == 'rot'))
# Deviance Null-Modell
-2 * log_likelihood
wein_glm
Modell bestimmen, nur musst du hier die eigentlichen vom Modell vorhergesagten Werte verwenden, die du mit fitted()
extrahieren kannst.# Wahrscheinlichkeit weiss für wein_glm
wsk_weiss <- fitted(wein_glm)
# Log likelihood wein_glm
log_likelihood <- sum(log(wsk_weiss[wein$Farbe == 'weiss'])) +
sum(log(1-wsk_weiss[wein$Farbe == 'rot']))
# Deviance wein_glm
-2 * log_likelihood
Die Deviance Werte bzw. die log-Likelihoods lassen sich nur schwer absolut beurteilen. Das bedeutet auch, dass nicht ohne weiteres klar ist, ob wein_glm
tatsächlich gut die Farbe
erklärt. Bei standard Regression ist das leichter, da R-squared
notwendigerweise zwischen 0
und 1
liegt. Um die Interpretation zu erleichtert, bietet R daher noch das Akaike Informationskriterium (AIC
), welches die Deviance in Vergleich zur Anzahl der Parameter setzt, d.h., AIC <- Deviance + 2*k
, wobei k
die Anzahl der Parameter (Regressionsgewichte ist). Wenn das AIC
niedriger ist als die Null-Deviance, dann kann man klar davon ausgehen, dass das Modell (wein_glm
) besser ist als das Null-Modell. Hier ist dies der Fall, was nicht überraschend ist, da die Prädiktoren jeweils sehr starke Effekte aufweisen.
Dass starke Effekte vorliegen, lässt sich über sogenannte Odds Ratios (OR) darstellen. ORs sind eine übliche Effektstärke in der logistischen Regression, die sich leicht aus Regressionsgewichten ableiten lässt und zwar gilt OR = exp(b)
, wobei b
das Regressionsgewicht ist. Berechne die ORs für die drei im Modell enthaltenen Prädiktoren. Es ist i.d.R. einfacher mit positiven Gewichten zu arbeiten, solange man die Richtung in Erinnerung behält.
# Odds ratios
exp(abs(XX))
exp(abs(XX))
exp(abs(XX))
p / (1-p)
) dafür, dass ein Wein rot ist, um das 2.59
fache steigen, wenn man Alkohol um 1 Volumenprozent zunimmt, um das 1.33e+279
(fast unzählbar) fache, wenn Dichte um einen Wert von 1 steigt, und um das 2477
fache, wenn die Sulphatkonzentration um 1 steigt. Ziemlich extreme Werte oder? Könnte vielleicht wiederum mit der Skalierung zu tun haben. Rechne die logistische Regression erneut, diesmal mit skalierten Variablen.# Skalierungsfunktion
skaliere = function(x) (x - mean(x))/sd(x)
# Logistische Regression mit Skalierung
wein_glm <- glm(formula = XX ~ XX + XX + XX,
data = XX %>% mutate_if(is.numeric, skaliere),
family = 'XX')
# Odds ratios
exp(abs(XX))
exp(abs(XX))
exp(abs(XX))
Dichte
war, dass die Werte in Dichte
nur zwischen 0.987
und 1.039
, was bedeutet, dass eine Veränderung von 1 eine Veränderung weit über den eigentlichen Wertebereich der Variablen hinausgeht. Wiederum, wenn es um den Vergleich der Gewichte geht und diese leicht interpretierbar werden sollen, erwäge zu standardisieren.bigmart
Datensatz. In diesem Abschnitt geht es Verkäufe
(Kriterium) auf Basis des maximalen Preis (Max_Preis
) und der Visibilität
vorherzusagen. Da es sich beim Kriterium um Häufigkeiten dreht empfiehlt sich eine Regression mit family = 'poisson'
. Verwende das Template um eine solche zu fitten.# Poisson Regression
bigmart_glm <- glm(formula = XX ~ XX + XX,
data = XX,
family = 'XX')
Printe zunächst das Modell. Wie interpretierst du den Output?
Der Output ist genauso aufgebaut wie bei der logistischen Regression. Entgegen was man intuitiv meinen könnte zeigt sich das Preis einen positiven Einfluss auf die Verkäufe hat und Visibilität einen negativen. Die Gewichte der Poisson-Regression sind ebenfalls per exp()
-Funktion zu transformieren, um sie richtig interpretieren zu können. In diesem Fall ist es aber nicht einfach die Gewichte separat von einander zu betrachten. Siehe die Beispiele unten.
# Interpretation der Gewichte
exp(6.71409) # Erwartete Häufigkeit wenn Max_Preis = 0 & Visibilität = 0
exp(6.71409 + 0.00715) # Erwartete Häufigkeit wenn Max_Preis = 1 & Visibilität = 0
exp(6.71409 - 2.11849) # Erwartete Häufigkeit wenn Max_Preis = 0 & Visibilität = 1
exp(6.71409 + 0.00715 - 2.11849) # Erwartete Häufigkeit wenn Max_Preis = 1 & Visibilität = 1
Evaluiere nun die summary()
. Was beobachtest du?
Alle Effekte sind signifikant, trotz des teilweise kleinen Einflusses (im Fall von Max_Preis
), was sicherlich der grossen Stichprobe zu schulden ist. Nichtsdestotrotz, zeigt die Residual deviance
einen klar niedrigeren Wert als die Null deviance
. Was sagt dir dies?
glm.nb()
aus dem MASS
Paket um eine solche Regression mit denselben Variablen wie zuletzt zu rechnen. Verwende das Template unten.# Lade MASS
library(MASS)
# Negative Binomial Regression
bigmart_glm.nb <- glm.nb(formula = XX ~ XX + XX,
data = XX)
Printe das Objekt und vergleiche die Regressionsgewichte mit denen der Poisson Regression. Hat sich was verändert?
Nein, die Regressionsgewichte sind sehr stabil geblieben. Print die summary()
. Was beobachtest du?
Alle Effekt sind immer noch signifikant. Was sich jedoch geändert hat sind die Deviances, was mit einem komplexen zweistufigen Fitting-Prozess zu tun hat. Vergleichbar sind aber immer noch AIC und der Wert der ganz am Ende angegeben wird, welcher die volle Deviance für das Modell ist. Du wirst feststellen, dass beide Werte deutlich besser ausfallen, als bei der Poisson Regression, was nahelegt, dass das Negative-Binomial Modell tatsächlich geeigneter ist, die Daten zu beschreiben.
# 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
Datei | Zeile | Spalte |
---|---|---|
wein.csv | 6497 | 13 |
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 Weins von 1-9 |
Farbe | Roter oder weisser Wien |
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 sauerer. |
Sulphate | Sulphatkontration im Wein |
Alkohol | Alkoholkonzentration im Wein in % |
Der bigmart.csv
Datensatz enthält Verkaufszahlen einer nationalen Supermarktkette in Nepal aus dem Jahre 2013 für 1559 Produkte und 10 verschiedene Supermarktfilialen.
Name | Beschreibung |
---|---|
Verkäufe | Verkaufszahlen |
Gewicht | Gewicht des Produkts |
Fettgehalt | Normal oder niedrig (low-fat) |
Visibilität | Verkaufsfläche für das Produkt |
Typ | Art des Produkts |
Max_Preis | Maximaler Preis |
Ladengrösse | Ladengrösse |
Package | Installation |
---|---|
tidyverse |
install.packages("tidyverse") |
MASS |
install.packages("MASS") |
Function | Package | Description |
---|---|---|
glm |
stats |
Fitte ein generalisiertes lineares Modell |
fitted |
stats |
Extrahiere vorhergesagte Werte |
residuals |
stats |
Extrahiere Residuen |
gml.nb |
MASS |
Fit a (generalized) linear model with negative binomial links |
Hilfreiche Vignette über Regressionen für Häufigkeiten mit R.
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.