Logistic Regression

(1)

Logistic Regression

Jonathan Harrington

Bitte lost.txt und lost2.txt laden – siehe Befehle in logistic.txt

Bitte lost.txt und lost2.txt laden – siehe

Befehle in logistic.txt

(2)

1. Logistic Regression: allgemeine Einführung

Dalgaard, P. (2002) Introductory Statistics with R.

Insbesondere Kap. 11

D. Cook, P. Dixon, W. M. Duckworth, M.S. Kaiser, K. Koehler, W.

Q. Meeker and W. R. Stephenson. Binary Response and Logistic Regression Analysis.

http://www.faculty.sbc.edu/bkirk/Biostatistics/course%20documents%20for

%202006/Logistic%20Regression%20Analysis.doc

Literatur

Baayen, R.H. Analyzing Linguistic Data: A practical introduction to Statistics. S. 213-234

Johnson, Keith (in press). Quantitative Methods in Linguistics.

Blackwell. Kapitel 5.

Verzani, J. (2005). Using R for Introductory Statistics

(Ebook ueber die LMU UB). Kapitel 12

(3)

1. Logistic Regression: allgemeine Einführung

Mit logistic Regression wird eine Regressionslinie an Proportionen angepasst.

Aus verschiedenen Gründen kann jedoch die lineare (least-squares) Regression nicht auf Proportionen angewandt werden.

Vor allem liegen Proportionen zwischen 0 und 1 während lineare Regression keine solchen Grenzen kennt (und

daher könnte ein lineares Regressionsmodell Proportionen unter 0 oder über 1 vorhersagen).

Außerdem wird in der linearen Regression eine konstante

Varianz angenommen; jedoch kann bewiesen werden,

dass je höher der Proportionsdurchschnitt, umso größer

die Varianz.

(4)

Diese (und andere) Probleme können überwunden werden:

1. wenn log-odds statt Proportionen modelliert werden

2. Durch Einsetzung von 'maximum likelihood' anstatt 'least squares'.

1. Logistic Regression: allgemeine Einführung

y = mx + b

logodds(y) = mx + b

Least-squares Regression Logistic Regression

Es wird nicht angenommen, dass die Werte Stichproben aus einer Normalverteilung sind.

Ein Vorteil von logistic Regression:

(5)

Einige Daten

high low 1950 30 5 1960 18 21 1971 15 26 1980 13 20 1993 4 32 2005 2 34

In 1950 produzierten 30 Sprecher /lo:st/ und 5 /lɔst/.

jahr = as.numeric(rownames(lost)) jahr = jahr - 1950

lost

(6)

Log-odds

p q=1-p Odds = p/q Log-Odds = log(p/q)

0.8 0.2 log(4) = 1.39

0.5 0.5 1 0

4 bedeutet 4:1 (wie im Pferderennen). Die Wahrscheinlichkeit vom Erfolg (p) ist 4 Mal so groß wie Scheitern (q)

p: Proportion 'Erfolg'.

lo:st lɔst

32 8

n 40

p 0.8 (prop. lOst)

(prop. lo:st)

(7)

Log-odds also log (p/q) als Funktion von p Log-odds haben Werte zwischen ±∞

Log-odds

(8)

2. Anwendung der logistic Regression in R: glm()

Das Ziel: nach der Anwendung von logistic Regression geben wir einen beliebigen Jahrgang ein, und das Modell soll uns die Proportion von /lo:st/ vorhersagen. zB Eingabe 1962, Proportion (lo:st) = ?

Jahr ist daher in diesem Fall die unabhängige Variable, Proportion von /lo:st/ die abhängige Variable.

Mit logistic Regression ist die abhängige Variable immer ein kategorialer Wert von 2 Möglichkeiten: ja oder nein, rot oder grün, 0 oder 1, weiblich oder

männlich, wach oder eingeschlafen, /lo:st/ oder /lɔst/,

Erfolg oder Scheitern, usw.

(9)

g = glm(lost ~ jahr, binomial) Ergebnis: ein Log-

Odd pro Jahr

unabhängige Variable (der Jahrgang)

bedeutet:

logistic

Regression ('binomial' weil wie in der

binomialen Verteilung wir mit 2 Werten (ja/nein,

Erfolg/Scheitern zu tun haben).

wird modelliert durch Abhängige Variable

Eine 2-spaltige Matrix: Anzahl von 'ja' und 'nein' (hier /lo:st/ und /lOst/)

lost

high low

1950 30 5

1960 18 21

1971 15 26

1980 13 20

1993 4 32

2005 2 34

(10)

Da die Ausgabe der Regression in log-odds ist,

müssen wir die Proportionen ebenfalls umwandeln, wenn wir die Regressionslinie sehen wollen.

Die vorhergesagten Werte überlagern

text(jahr, predict(g), "x", col=3)

Eine Abbildung der Daten in diesem Raum:

# Proportion von /lo:st/ berechnen p = lost[,1]/apply(lost, 1, sum)

# log-odds

lodd = log(p/(1-p))

plot(jahr, lodd, type="b")

# Regressionslinie überlagern abline(g, col=2)

3. Abbildung der Regressionslinie

(11)

Vorhersage: Wert für 1962

abline(h=m * 12 + k, lty=2, col="blue") abline(v= 12)

m * x + k m * x + k coef(g)

(Intercept) jahr 1.10432397 -0.07026313

k = coef(g)[1]

m = coef(g)[2]

(12)

Wir können durch die Transformation (2) die Regressionslinie auch in einem Raum von Jahr x Proportionen abbilden.

Von Proportionen in log-odds

Von log-odds zurück in Proportionen

1 )

log( p L p

  L = log(p/(1-p))

L L

e p e

 

1 p = exp(L)/(1+exp(L)) [1] 1.386294

[1] 0.8

Abbildung der Regression

(1)

(2)

p = 0.8

(13)

curve(exp(mx + k)/(1+ exp(mx+k)), xlim=c(0, 60), add=T, col=2)

# Regression überlagern

# Proportionen von /lo:st/ berechnen p = lost[,1]/apply(lost, 1, sum)

plot(jahr,p)

# Abbildung Jahr x Proportionen

Die Regression

) (

1 ^mx ^k

k mx

e

y e _



 

coef(g)

(Intercept) jahr 1.10432397 -0.07026313

k = coef(g)[1]

m = coef(g)[2]

Abbildung: Jahr x Proportionen

(14)

die vorhergesagten Werte liegen auf der Kurve:

vorher = predict(g)

text(jahr, exp(vorher)/(1+exp(vorher)), "x", col=3)

Abbildung Jahr x Proportionen

1962 exp(m12 + k)/(1+ exp(m12+k))

) (

1 ^mx ^k

k mx

e

y e _



 

(15)

Lineare Regression: R

²

oder adjusted R

²

und ein F-test Logistic Regression: G

²

und ein 

²

-test.

G

²

= Null deviance – residual deviance

Für ein signifikantes Ergebnis wollen wir daher, dass Null deviance hoch und Residual deviance klein ist.

wenn dieser Wert 0 wäre, dann wären alle Proportionen in allen Jahren gleich (und die Regressionslinie wäre horizontal)

je höher dieser Wert, umso unwahrscheinlicher ist es, dass die Werte überhaupt durch die Regression

modelliert werden können.

Signifikanz-Test

Was ist die Wahrscheinlichkeit, dass die Proportion von /lo:st-

lɔst/ durch den Jahrgang vorhergesagt werden kann?

(16)

69.3634 - 8.2422

[1] 61.1212

G

²

= Null deviance - residual deviance summary(g)

Null deviance: 69.3634 on 5 degrees of freedom

Residual deviance: 8.2422 on 4 degrees of freedom

…

Der Test mit anova() ist ob G

²

signifikant von 0 abweicht:

anova(g, test="Chisq")

Df Deviance Resid. Df Resid. Dev P(>|Chi|) NULL 5 69.363 jahr 1 61.121 4 8.242 5.367e-15

Eine Analyse der logistischen Regression mit Proportion als

abhängige Variable und Jahr als unabhängige Variable zeigte, dass die Proportionen einem Trend folgen (

²

(1)=61.2, p < 0.001)

g = glm(lost ~ jahr, binomial)

(17)

logistische Regression logistische Regression

dd = read.delim(paste(pfad, "kj.txt", sep="/"))

lneu = read.table(paste(pfad, "lneu.txt", sep="/"))

Bitte kj.txt und lneu.txt herunterladen und diese

Befehle aus logbefehle2.txt eingeben

(18)

g = glm(var ~ Faktor, binomial) Ergebnis: ein Log-Odd pro Jahr

Abhängige Variable Eine 2-spaltige

Matrix Anzahl von 'ja' oder 'nein'

high low 1950 30 5 1960 18 21 1971 15 26 1980 13 20 1993 4 32 2005 2 34

Oder ein Vektor

high high low high high high low ...

Faktor

Trend Analyse

nimmt der high/low Anteil mit

zunehmenden Jahrgang ab?

Kategoriale Analyse Unterscheiden sich

männlich und weiblich in der Anteil von high/low?

Logistische Regression

(19)

Kategoriale Analyse (Ein Faktor) Kategoriale Analyse (Ein Faktor) dd = read.delim(paste(pfad, "kj.txt", sep="/")) head(dd)

names(dd) 1. Tabelle 2. Barplot

3. Modell: logistische Regression 4. Signifikanz

tab = with(dd, table(str, age))

Frage: Wird die Verwendung von /str/ vs. /Str/ vom Alter beeinflusst?

barplot(tab, beside=T, legend=T)

anova(g, test="Chisq")

g = glm(str ~ age, binomial, data=dd) Kategorial

5. Ergebnis

Die str/Str Verteilung wird signifikant vom Alter beinflusst

(

^

[2] = 17.4, p < 0.001)

(20)

Trend Analyse (ein numerischer Faktor) Trend Analyse (ein numerischer Faktor) lneu = read.table(paste(pfad, "lneu.txt", sep="/")) head(lneu)

names(lneu) 1. Tabelle

2. Jahr x Proportion Abbildung

3. Modell: logistische Regression 4. Signifikanz

Frage: Wird die Verwendung von high vs. low vom Jahrgang beeinflusst?

5. Ergebnis

p = tab[1,]/apply(tab, 1, sum) zeit = as.numeric(colnames(tab)) tab = with(lneu, table(Vowel, Jahr))

g = glm(Vowel ~ Jahr, binomial, data=lneu) anova(g, test="Chisq")

Numerisch

Die high/low Verteilung wird signifikant vom Jahrgang beinflusst (

^

[1] = 61.1, p < 0.001)

plot(zeit, p, type="b", xlab="Jahr (Re. 1950)", ylab="Proportion HIGH")

(21)

Zwei Faktoren.

high low 00.0 16 0 10.0 9 6 21.0 8 10 30.0 8 7 43.0 4 10 55.0 1 15 00.1 14 5 10.1 9 15 21.1 7 16 30.1 5 13 43.1 0 22 55.1 1 19

M

W

1950 1960 1971 1980 1993 2005 1950 1960 1971 1980 1993 2005

} }

(a) Gibt es einen

Trend? Also weniger high in späteren

Jahren?

(b) Ist die Proportion high/low in M und W unterschiedlich

verteilt?

Interaktion: Jahr x Geschlecht

Faktor: Jahr

Faktor: Geschlecht

Ist der Trend derselbe für M und W?

beide = paste(G, Jahr, sep=".")

tab = table(beide, Vowel)

(22)

Zuerst eine Abbildung…

temp = Vowel == "high"

interaction.plot(Jahr, G, temp)

Nimmt die Proportion von high in späteren Jahren ab? (Die Unterschiede zwischen m und w ignorieren).

Ja Nein Vielleicht

Unterscheiden sich m und w in der Proportion von high?

(Die Unterschiede in den Jahrgängen ignorieren).

Ja Nein Vielleicht

attach(lneu)

(23)

Coefficients:

(Intercept) Jahr G

1.87754 -0.07524 1.20282

logodds(lo:st) = 1.87754 - 0.07524 * Jahr + 1.20282 * G mehrg = glm(Vowel ~ Jahr + G, binomial, data=lneu)

anova(mehrg, test="Chisq")

Df Deviance Resid.Df Resid. Dev P(>|Chi|) NULL 219 290.572

J 1 61.121 218 229.451 5.367e-15 G 1 12.822 217 216.629 3.425e-04

M und W unterscheiden sich in der Proportion von high/low, 

²

(1) = 12.82, p < 0.001

Die Proportion von high nimmt in späteren Jahren ab,



²

(1) = 61.12, p < 0.001.

mehrg

(24)

Mit 2 oder mehr Variablen soll auch geprüft werden, ob sie miteinander interagieren.

Eine Interaktion zwischen den unabhängigen Variablen – in diesem Fall Geschlecht und Jahrgang – liegt vor, wenn sie eine unterschiedliche Wirkung auf die abhängige

Variable ausüben wie in 1 und 2, aber nicht in 3 und 4

Die Interaktion zwischen 2 Variablen

1950 2000

prop(high)

1950 2000 1950 2000 1950

2000

prop(high) prop(high) prop(high)

1 2 3 4

m

w

(25)

Wenn eine Interaktion vorliegt, dann können signifikante Ergebnisse in einer der unabhängigen Variablen nicht uneingeschränkt akzeptiert werden.

zB wenn eine Interaktion vorkommt, gibt es vielleicht eine

Wirkung von Jahrgang auf die Proportion von high/low nur in Männern aber nicht in Frauen usw.

Die Interaktion zwischen 2 Variablen

dies scheint aber hier nicht

der Fall zu sein.

(26)

Die Interaktion zwischen 2 unabhängigen Variablen, A und B, kann in R mit A:B geprüft werden.

Daher in diesem Fall

g = glm(Vowel ~ Jahr + G + Jahr:G, binomial, data=lneu) Eine Abkürzung dafür (und mit genau demselben Ergebnis)

g = glm(Vowel ~ Jahr * G, binomial, data=lneu)

Die Interaktion zwischen 2 Variablen

anova(g, test="Chisq")

Df Deviance Resid. Df Resid. Dev P(>|Chi|) NULL 219 290.572 Jahr 1 61.121 218 229.451 5.367e-15 ***

G 1 12.822 217 216.629 0.0003425 ***

Jahr:G 1 0.017 216 216.612 0.8962973

d.h. die Interaktion ist nicht signifikant und J:G kann

aus dem Regressionsmodell weggelassen werden.

(27)

library(MASS) stepAIC(g)

Dies wird auch durch stepAIC() bestätigt:

AIC wird kleiner wenn wir Jahr:G weglassen

Start: AIC=224.61 Vowel ~ Jahr * G

Df Deviance AIC - Jahr:G 1 216.63 222.63

<none> 216.61 224.61

daher bleiben wir bei:

Call: glm(formula = Vowel ~ Jahr + G, family = binomial, data = lneu)

(28)

pfad = "Pfad wo Sie kj.txt gespeichert haben"

kj = read.delim(paste(pfad, "kj.txt", sep="/")) Inwiefern wird die Wahl zwischen 'str' und 'Str' durch

Alter und Herkunft (age und region) Alter und soziale Klasse (age, class) Alter und Geschlecht (age, gender)

beeinflusst?