Logistic Regression

1. Logistic Regression Jonathan Harrington Befehle: logistic.txt

2. 1. Logistic Regression: allgemeine Einführung Literatur Baayen, R.H. Analyzing Linguistic Data: A practical introduction to Statistics. S. 213-234 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 Dalgaard, P. (2002) Introductory Statistics with R. Insbesondere Kap. 11 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 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. Mit logistic Regression wird eine Regressionslinie an Proportionen angepasst.

4. 1. Logistic Regression: allgemeine Einführung Logistic Regression logodds(y) = mx + b Least-squares Regression y = mx + b Ein Vorteil von logistic Regression: Es wird nicht angenommen, dass die Werte Stichproben aus einer Normalverteilung sind. 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'.

5. Einige Daten lost In 1950 produzierten 30 Sprecher /lo:st/ und 5 /lɔst/. high low 1950 30 5 1960 18 21 1971 15 26 1980 13 20 1993 4 32 2005 2 34 jahr = as.numeric(rownames(lost)) jahr = jahr - 1950

6. Log-odds 4 bedeutet 4:1 (wie im Pferderennen). Die Wahrscheinlichkeit vom Erfolg (p) ist 4 Mal so groß wie Scheitern (q) 0.5 0.5 1 0 p: Proportion 'Erfolg'. p lo:st lɔst n 32 8 40 0.8 (prop. lo:st) (prop. lOst) p q=1-p Odds = p/q Log-Odds = log(p/q) 0.8 0.2 log(4) = 1.39

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

8. 2. Anwendung der logistic Regression in R: glm() 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. 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.

9. Ergebnis: ein Log-Odd pro Jahr unabhängige Variable (der Jahrgang) wird modelliert durch bedeutet: logistic Regression ('binomial' weil wie in der binomialen Verteilung wir mit 2 Werten (ja/nein, Erfolg/Scheitern zu tun haben). 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 g = glm(lost ~ jahr, binomial)

10. 3. Abbildung der Regressionslinie 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) Die vorhergesagten Werte überlagern text(jahr, predict(g), "x", col=3) Da die Ausgabe der Regression in log-odds ist, müssen wir die Proportionen ebenfalls umwandeln, wenn wir die Regressionslinie sehen wollen.

11. Vorhersage: Wert für 1962 neuerwert = data.frame(jahr = 12) ergebnis = predict(g, neuerwert, se.fit=T) ergebnis$fit abline(h=ergebnis$fit, lty=2, col="blue")

12. Abbildung der Regression Wir können durch die Transformation (2) die Regressionslinie auch in einem Raum von Jahr x Proportionen abbilden. Von Proportionen in log-odds p = 0.8 (1) L = log(p/(1-p)) [1] 1.386294 Von log-odds zurück in Proportionen (2) p = exp(L)/(1+exp(L)) [1] 0.8

13. coef(g) (Intercept) jahr 1.10432397 -0.07026313 m = coef(g)[2] k = coef(g)[1] Abbildung: Jahr x Proportionen # Proportionen von /lo:st/ berechnen p = lost[,1]/apply(lost, 1, sum) # Abbildung Jahr x Proportionen plot(jahr,p) Die Regression # Regression überlagern curve(exp(m*x + k)/(1+ exp(m*x+k)), xlim=c(0, 60), add=T, col=2)

14. Abbildung Jahr x Proportionen und die vorhergesagten Werte liegen wieder auf der Linie: vorher = predict(g) text(jahr, exp(vorher)/(1+exp(vorher)), "x", col=3)

15. 4. Anwendung der logistic Regression auf die Berechung einer perzeptiven Grenzen zwischen Kategorien Experiment. Anhand der Sprachsynthese wurde ein F2-Kontinuum in 11 Schritten synthetisiert. 5 Vpn. (L1-Englisch) mussten zu jedem Stimulus mit "I" oder "U" antworten. Bei welchem F2-Wert liegt die Grenze zwischen den Vokalen? Die Anzahle der Bewertungen ist hier: ui u i 2311 0 5 2176 0 5 2023 0 5 1885 0 5 1770 0 5 1667 0 5 1548 2 3 1437 4 1 1351 5 0 1269 5 0 Bei F2 = 1437 Hz gab es 4 Urteil für "U", ein Urteil für "I"

16. Bei 1437 Hz waren 80% der Urteile "U" (und daher 20% "I") Ein Vektor von Proportionen p = ui[,1]/apply(ui, 1, sum) p 2311 2176 2023 1885 1770 1667 1548 1437 1351 1269 0.0 0.0 0.0 0.0 0.0 0.0 0.4 0.8 1.0 1.0 Eine Abbildung von F2 als Funktion dieser Proportionen f2werte = as.numeric(rownames(ui)) plot(f2werte, p, ylab="Proportion /u/ Urteile", xlab="F2 (Hz)")

17. Eine logistische Regression an diese Werte anpassen Die Urteile aus den F2-Werten vorhersagen logui = glm(ui ~ f2werte, family=binomial) Mit der Methode auf S. 13-14 die logistische Regressionskurve überlagern Die Koeffiziente m = coef(logui)[2] k = coef(logui)[1] Die logistische Regressionskurve curve(exp(m*x + k)/(1+ exp(m*x+k)), xlim=c(1200, 2400), add=T, col=2)

18. Die 50% Grenze (Umkipppunkt) = zu welchem F2-Wert, ist ein Urteil für "I" genauso wahrscheinlich wie ein Urteil für "U"? Es kann bewiesen werden, dass dies mit -k/m in dieser Formel gegeben wird (in diesem Beispiel ist y die Proportion, p, und x ist F2werte) abline(v=-k/m, lty=2, col="blue") -k/m 1516.723

20. 5. Signifikanz-Test je höher dieser Wert, umso unwahrscheinlicher ist es, dass die Werte überhaupt durch die Regression modelliert werden können. wenn dieser Wert 0 wäre, dann wären alle Proportionen in allen Jahren gleich (und die Regressionslinie wäre horizontal) Was ist die Wahrscheinlichkeit, dass die Proportion von /lo:st-lɔst/ durch den Jahrgang vorhergesagt werden kann? Lineare Regression: R2 oder adjusted R2 und ein F-test Logistic Regression: G2 und ein c2-test. G2 = Null deviance – residual deviance Für ein signifikantes Ergebnis wollen wir daher, dass Null deviance hoch und Residual deviance klein ist.

21. Der Test mit anova() ist ob G2 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 G2 = Null deviance - residual deviance g = glm(lost ~ jahr, binomial) summary(g) … Null deviance: 69.3634 on 5 degrees of freedom Residual deviance: 8.2422 on 4 degrees of freedom 69.3634 - 8.2422 [1] 61.1212 Die Proportionen folgen einem Trend (c2(1)=61.2, p < 0.001)

22. 6. Zwei unabhängige Variablen. von diesen 26 waren 10 von Männern und 16 von Frauen erzeugt. (a) Gibt es einen Trend? Also weniger [lo:st] in späteren Jahren? (b) Ist die Proportion [lost]/[lo:st] in M und F unterschiedlich verteilt? Hier sind genau dieselben Daten aber zusätzlich nach männlich-weiblich aufgeteilt. female lost 1950 1960 1971 1980 1993 2005 n n 16 9 8 8 4 1 y n 14 9 7 5 0 1 n y 0 6 10 7 10 15 y y 5 15 16 13 22 19 In 1971 waren 26 Tokens [lost] und 15 [lo:st] lost high low 1950 30 5 1960 18 21 1971 15 26 1980 13 20 1993 4 32 2005 2 34 8M, 7F

23. female lost 1950 1960 1971 1980 1993 2005 n n 16 9 8 8 4 1 y n 14 9 7 5 0 1 n y 0 6 10 7 10 15 y y 5 15 16 13 22 19 Dies ist ein Problem der mehrfachen Logistic Regression: (also in diesem Fall eine Linie im 3D-Raum) logodds (lo:st) = b0 + b1year + b2Geschlecht logodds(lo:st) (b0 ist das Intercept, b1 und b2 die Neigungen) Geschlecht Und eine gerade Linie in einem 3D-Raum Year

24. Daten-Vorbereitung J = c(jahr, jahr) J [1] 0 10 21 30 43 55 0 10 21 30 43 55 G = c(rep(0, 6), rep(1, 6)) G [1] 0 0 0 0 0 0 1 1 1 1 1 1 pfad = "das Verzeichnis wo ich lost2.txt gespeichert habe" lost2 = as.matrix(read.table(paste(pfad, "lost2.txt", sep="/"))) } high low 0.0 16 0 10.0 9 6 21.0 8 10 30.0 8 7 43.0 4 10 55.0 1 15 0.1 14 5 10.1 9 15 21.1 7 16 30.1 5 13 43.1 0 22 55.1 1 19 1950 1960 1971 1980 1993 2005 high = Spalte 1 = /lo:st/ M low = Spalte 2 = /lOst/ } 1950 1960 1971 1980 1993 2005 W

25. Zuerst eine Abbildung… Nimmt die Proportion von /lo:st/ in späteren Jahren ab? (Die Unterschiede zwischen m und f ignorieren). G m 0.8 Ja Nein Vielleicht f mean of p 0.4 Unterscheiden sich m und f in der Proportion von /lo:st/? (Die Unterschiede in den Jahrgängen ignorieren). 0.0 J 0 10 21 30 43 55 Ja Nein Vielleicht p = lost2[,1]/apply(lost2, 1, sum) interaction.plot(J, G, p)

26. Wenn wir übrigens G weglassen, dann müssten wir trotz der anderen Aufteilung der Daten das gleiche Ergebnis wir vorhin bekommen: g2 = glm(lost2 ~ J, binomial) anova(g2, test="Chisq") Analysis of Deviance Table Df Deviance Resid. Df Resid. Dev P(>|Chi|) NULL 11 89.557 year 1 61.121 10 28.436 5.367e-15 Modell berechnen… mehrg = glm(lost2 ~ J + G, binomial)

27. mehrg = glm(lost2 ~ J + G, binomial) anova(mehrg, test="Chisq") Df Deviance Resid.Df Resid. Dev P(>|Chi|) NULL 11 89.557 J 1 61.121 10 28.436 5.367e-15 G 1 12.822 9 15.613 3.425e-04 Die Proportion von 'lo:st' nimmt in späteren Jahren ab, c2(1) = 61.12, p < 0.001. M und F unterscheiden sich in der Proportion von lo:st/lOst, c2(1) = 12.82, p < 0.001 mehrg Coefficients: (Intercept) J Gm 1.87754 -0.07524 1.20282 Degrees of Freedom: 11 Total (i.e. Null); 9 Residual Null Deviance: 89.56 Residual Deviance: 15.61 AIC: 51.51 logodds(lo:st) = 1.87754 - 0.07524J+ 1.20282G

28. 7. Die Interaktion zwischen 2 Variablen prop(lo:st) prop(lo:st) prop(lo:st) prop(lo:st) f m 2000 1950 1950 2000 1950 2000 1950 2000 1 2 3 4 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

29. Die Interaktion zwischen 2 Variablen dies scheint aber hier nicht der Fall zu sein. 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 /lo:st/ nur in Männern aber nicht in Frauen usw.

30. Die Interaktion zwischen 2 Variablen Daher in diesem Fall g = glm(lost2 ~ J + G + J:G, binomial) Eine Abkürzung dafür (und mit genau demselben Ergebnis) g = glm(lost2 ~ J * G, binomial) anova(g, test="Chisq") Df Deviance Resid. Df Resid. Dev P(>|Chi|) NULL 11 89.557 J 1 61.121 10 28.436 5.367e-15 G 1 12.822 9 15.613 3.425e-04 J:G 1 0.017 8 15.596 0.896 d.h. die Interaktion ist nicht signifikant und J:G kann aus dem Regressionsmodell weggelassen werden. Die Interaktion zwischen 2 unabhängigen Variablen, A und B, kann in R mit A:B geprüft werden.

31. Dies wird auch durch stepAIC() bestätigt: library(MASS) stepAIC(g) AIC wird kleiner wenn wir J:G weglassen Start: AIC= 53.49 lost2 ~ J * G Df Deviance AIC - J:G 1 15.613 51.506 <none> 15.596 53.489 Df Deviance AIC <none> 15.613 51.506 - G 1 28.436 62.328 - J 1 80.018 113.910 Wir bleiben also bei Call: glm(formula = lost2 ~ J + G, family = binomial) Residual Deviance: 15.61 AIC: 51.51

32. Weitere Folien zum Durchlesen...

33. gmf = glm(lost2 ~ G, "binomial") anova(gmf, test="Chisq") Df Deviance Resid. Df Resid. Dev P(>|Chi|) NULL 11 89.557 G 1 9.539 10 80.018 0.002 M und F unterscheiden sich in der Proportion von lo:st/lOst (c2(1) = 9.5, p < 0.002). 8. Logistic Regression und zwei Ebenen Aus dem vorigen Beispiel wird auch klar, dass ähnlich wie c2 Logistic Regression angewandt werden kann, auch wenn die Gruppe nur aus 2 Ebenen besteht. Gibt es einen signifikanten Unterschied zwischen M und F?

34. high low 0 46 48 1 36 90 (kodiert nur nach M und F) = l.mf = c(0,1) gmf2 = glm(mf ~ l.mf, "binomial") anova(gmf2, test="Chisq") Wir bekommen dasselbe Ergebnis wenn Logistic Regression auf die entsprechende Tabelle angewandt wird: m = apply(lost2[1:6,], 2, sum) f = apply(lost2[7:12,], 2, sum) mf = rbind(m, f) rownames(mf) = c(0, 1) colnames(mf) = c("high", "low") mf lost2 high low 0.0 16 0 10.0 9 6 21.0 8 10 30.0 8 7 43.0 4 10 55.0 1 15 0.1 14 5 10.1 9 15 21.1 7 16 30.1 5 13 43.1 0 22 55.1 1 19

35. und man bekommt dann fast das gleiche Ergebnis mit einem c2-Test, der direkt auf die Tabelle angewandt wird: chisq.test(mf) Pearson's Chi-squared test with Yates' continuity correction data: mf X-squared = 8.6985, df = 1, p-value = 0.003185 Ein c2-Test kann jedoch nicht verwendet werden, bei einer Gruppenanzahl von > 2 …

36. lost3 high low alt.0 43 35 alt.1 30 15 jung.0 3 13 jung.1 6 75 Hier sind dieselben Daten aufgeteilt in 2 Altersgruppen sowie M/F Gruppe 1 = Vokal = high/low Gruppe 2 = Geschl = M/F (=0/1) Gruppe 3 = Alter = alt/jung 3 Gruppen jeweils 2 Ebenen Haben (a) Alter und (b) Geschlecht einen Einfluss auf die Proportion von /lo:st/?

37. Zuerst eine Abbildung high low alt.0 43 35 alt.1 30 15 jung.0 3 13 jung.1 6 75 # Alter kodieren A = c(0, 0, 1, 1) # Geschlecht kodieren G = c(0, 1, 0, 1) prop = lost3[,1]/apply(lost3, 1, sum) interaction.plot(A, G, prop)

38. G 0.6 0 1 0.5 0.4 mean of prop 0.3 0.2 0.1 A 1 0 Signifikanter Einfluss auf lo:st/lOst? ja nein vielleicht im Alter? im Geschlecht? ja nein vielleicht ja nein vielleicht Interaktion zwischen A und G?

39. g = glm(lost3 ~ A * G, binomial) anova(g, test="Chisq") Df Deviance Resid. Df Resid. Dev P(>|Chi|) NULL 3 67.758 A 1 64.452 2 3.306 9.893e-16 G 1 0.398 1 2.908 0.528 A:G 1 2.908 0 -4.441e-16 0.088 Es gab einen signifikanten Einfluss vom Alter (c2(1)=64.2, p < 0.001) aber nicht vom Geschlecht auf die Proportion von /lo:st/. Die Interaktion zwichen Alter und Geschlecht war nicht signifikant (p > 0.05).