Logistische Regression - Beispiel in R

Der erste Teil der Artikelserie zur logistischen Regression stellt die logistische Regression als Verfahren zur Modellierung binärer abhängiger Variablen vor. Der zweite Teil geht auf Methoden für die Beurteilung der Klassifikationsgüte ein. In diesem Artikel wird nun die Anwendung des Verfahrens an einem konkreten Beispiel, der Klassifikation von Weinen, mithilfe der Statistik-Software R gezeigt.

Datensatz: Klassifikation von Weinen

Es beschäftigt uns (wie schon bei der Vorstellung der linearen Regression) auch bei der logistischen Regression wieder das Thema Wein. Diesmal geht es jedoch nicht darum, die Qualität des Weines mittels Regression zu bestimmen, vielmehr soll nun anhand der chemischen Eigenschaften des Weins seine Farbe (rot oder weiß) bestimmt werden. Der verwendete Datensatz enthält insgesamt 6497 Beobachtungen, davon gehören 1599 zu den Rot- und 4898 zu den Weißweinen. In der Spalte "color" wird die Farbe spezifiziert, wobei 0 für "rot" und 1 für "weiß" steht. Die Variable "quality" enthält eine Einschätzung der Qualität des Weines auf einer Skala von 0 bis 10. Darüber hinaus geben 11 weitere Variablen Aufschluss über die chemischen Eigenschaften der Weine:

• color (0=rot, 1=weiß), quality (zwischen 0 und 10), fixed acidity, volatile acidity, citric acid, residual sugar, chlorides, free sulfur dioxide, total sulfur dioxide, density, pH, sulphates, alcohol

 # Rotweindatensatz einlesen 
red <- read.csv2("http://archive.ics.uci.edu/ml/machine-learning-databases/wine-quality/winequality-red.csv", dec = ".", header = TRUE) 

# Weißweindatensatz einlesen
white <- read.csv2("http://archive.ics.uci.edu/ml/machine-learning-databases/wine-quality/winequality-white.csv", dec = ".", header = TRUE) 

# jedem der beiden Datensätze eine Spalte "color" mit 0 bei Rotweinen und 1 bei Weißweinen anfügen
red$color <- 0
white$color <- 1 

# Zusammenführen der zwei Datensätze zu einem Datensatz "wine"
wine <- rbind(red, white) 

# Löschen der nun überflüssigen Einzeldatensätze
rm(list = c("red", "white"))

Modellierung mittels Logit

Im ersten Schritt verschaffen wir uns einen Überblick über den Datensatz und schätzen dann ein Logit-Modell mit allen zur Verfügung stehenden Variablen. Außer fixed.acidity und pH sind alle Variablen zu einem Niveau von signifikant. Als Beispiel für eine Interpretation wird der Regressionskoeffizient der Variable für den Gehalt der Zitronensäure citric.acid herangezogen. Steigt der Zitronensäuregehalt um eine Einheit (in diesem Fall 0.01), so erhöht sich die Chance

um den Faktor .

Dieses erste Modell wird im Folgenden auf seine Klassifikationsgüte überprüft.

Um herauszufinden ob die insignifikanten Variablen aus dem Modell entfernt werden können, bietet sich die Verwendung der step()-Funktion an (siehe R Code). Diese minimiert schrittweise (durch hinzufügen oder weglassen von Variablen) das Akaike Informationskriterium (AIC) des Modells. Tatsächlich bieten die Variablen fixed.acidity und pH keinen zusätzlichen signifikanten Erklärungsgehalt für unser Modell und können damit von der Analyse ausgeschlossen werden.

# einen Überblick über den Datensatz verschaffen: 
str(wine)

'data.frame':   6497 obs. of  13 variables:
 $ fixed.acidity       : num  7.4 7.8 7.8 11.2 7.4 7.4 7.9 7.3 7.8 7.5 ...
 $ volatile.acidity    : num  0.7 0.88 0.76 0.28 0.7 0.66 0.6 0.65 0.58 0.5 ...
 $ citric.acid         : num  0 0 0.04 0.56 0 0 0.06 0 0.02 0.36 ...
 $ residual.sugar      : num  1.9 2.6 2.3 1.9 1.9 1.8 1.6 1.2 2 6.1 ...
 $ chlorides           : num  0.076 0.098 0.092 0.075 0.076 0.075 0.069 0.065 0.073 0.071 ...
 $ free.sulfur.dioxide : num  11 25 15 17 11 13 15 15 9 17 ...
 $ total.sulfur.dioxide: num  34 67 54 60 34 40 59 21 18 102 ...
 $ density             : num  0.998 0.997 0.997 0.998 0.998 ...
 $ pH                  : num  3.51 3.2 3.26 3.16 3.51 3.51 3.3 3.39 3.36 3.35 ...
 $ sulphates           : num  0.56 0.68 0.65 0.58 0.56 0.56 0.46 0.47 0.57 0.8 ...
 $ alcohol             : num  9.4 9.8 9.8 9.8 9.4 9.4 9.4 10 9.5 10.5 ...
 $ quality             : int  5 5 5 6 5 5 5 7 7 5 ...
 $ color               : num  0 0 0 0 0 0 0 0 0 0 ...

# Definition des Modells 
modell <- as.formula("color ~ fixed.acidity + volatile.acidity + 
                      citric.acid + residual.sugar + chlorides + 
                      free.sulfur.dioxide + total.sulfur.dioxide + 
                      density + pH + sulphates + alcohol + quality")

logit <- glm(modell, family = binomial, data = wine) # Ausgabe der Schätzergebnisse 
summary(logit)

Call:
glm(formula = modell, family = binomial, data = wine)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-5.6178   0.0012   0.0188   0.0582   6.8678  

Coefficients:
                        Estimate  Std. Error z value             Pr(>|z|)    
(Intercept)           1875.95755   186.79663  10.043 < 0.0000000000000002 ***
fixed.acidity            0.40048     0.23335   1.716               0.0861 .  
volatile.acidity        -6.72197     1.06072  -6.337       0.000000000234 ***
citric.acid              2.61720     1.18468   2.209               0.0272 *  
residual.sugar           0.95622     0.10119   9.449 < 0.0000000000000002 ***
chlorides              -22.01150     3.98449  -5.524       0.000000033082 ***
free.sulfur.dioxide     -0.06080     0.01456  -4.177       0.000029600549 ***
total.sulfur.dioxide     0.05229     0.00499  10.479 < 0.0000000000000002 ***
density              -1875.04266   190.43498  -9.846 < 0.0000000000000002 ***
pH                       1.95933     1.42434   1.376               0.1689    
sulphates               -2.69271     1.24900  -2.156               0.0311 *  
alcohol                 -1.79214     0.27949  -6.412       0.000000000143 ***
quality                 -0.43387     0.20411  -2.126               0.0335 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 7250.98  on 6496  degrees of freedom
Residual deviance:  424.23  on 6484  degrees of freedom
AIC: 450.23

Number of Fisher Scoring iterations: 9

# Minimierung des AIC mittels step()-Funktion 
logitMinAIC <- step(logit) 
summary(logitMinAIC)

Call:
glm(formula = color ~ volatile.acidity + citric.acid + residual.sugar + 
    chlorides + free.sulfur.dioxide + total.sulfur.dioxide + 
    density + sulphates + alcohol + quality, family = binomial, 
    data = wine)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-5.6466   0.0010   0.0173   0.0553   6.1056  

Coefficients:
                         Estimate   Std. Error z value             Pr(>|z|)    
(Intercept)           1645.421879   121.145643  13.582 < 0.0000000000000002 ***
volatile.acidity        -7.102440     1.029847  -6.897     0.00000000000533 ***
citric.acid              2.831430     1.089867   2.598              0.00938 ** 
residual.sugar           0.871605     0.086265  10.104 < 0.0000000000000002 ***
chlorides              -24.384586     3.772300  -6.464     0.00000000010189 ***
free.sulfur.dioxide     -0.058600     0.014598  -4.014     0.00005965241437 ***
total.sulfur.dioxide     0.052241     0.004991  10.467 < 0.0000000000000002 ***
density              -1635.753541   120.408807 -13.585 < 0.0000000000000002 ***
sulphates               -3.056311     1.193103  -2.562              0.01042 *  
alcohol                 -1.560248     0.229300  -6.804     0.00000000001015 ***
quality                 -0.410699     0.199857  -2.055              0.03988 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 7250.98  on 6496  degrees of freedom
Residual deviance:  427.23  on 6486  degrees of freedom
AIC: 449.23

Number of Fisher Scoring iterations: 9

Beurteilung der Klassifikationsgüte im Logit

Zuerst wird eine Klassifikationstabelle erstellt, um zu erkennen wie viele Weine das Modell mit einem Schwellenwert von 0.5 (Standard) der richtigen Farbe zuordnet:

Es ist zu erkennen, dass 1580 der 1599 Rotweine und 4887 der 4898 Weißweine korrekt klassifiziert werden. Die Korrektklassifikationsrate beträgt . Das Modell kann damit als sehr gut angesehen werden. In der unteren Grafik ist darüber hinaus die ROC-Kurve für das Logit-Modell abgetragen. Die Fläche zwischen der Diagonalen und der ROC Kurve hat annähernd die maximale Größe. Eine weitere Anpassung des Modells an den Datensatz würde im Folgenden kaum eine Möglichkeit zur Optimierung des Klassifikationsergebnisses bieten.

Auch die Pseudo-Bestimmtheitsmaße deuten in diesem Fall auf eine gute Modellanpassung hin: Das McFadden R² beträgt 0.94, wobei bereits Werte ab 0.2/0.3 auf ein gutes Modell schließen lassen. Das Pseudo-Bestimmtheitsmaß von Nagelkerke hat hingegen den Vorteil, dass es analog zum Bestimmtheitsmaß bei der linearen Regression interpetiert werden kann. Nach Nagelkerke ergibt sich ein Wert von 0.97. Zuletzt bleibt noch das Cox&Snell R², dieses beträgt für unser Wein-Modell 0.65.

# Die benötigten Pakete laden
library(caret)
library(ggplot2)
library(plotROC)

# Erstellung einer Klassifikationstabelle
pred <- ifelse(fitted(logit) > 0.5, 1, 0)
confusionMatrix(factor(wine$color), factor(pred))

# ROC Kurve
df_predictions <- data.frame(color = wine$color, pred = fitted(logit))

ggplot(df_predictions, aes(d = color, m = pred)) + 
  geom_roc(n.cuts = 0) +
  geom_abline(slope = 1, intercept = 0, lty = 2)