close

Anmelden

Neues Passwort anfordern?

Anmeldung mit OpenID

Einführung in die Statistik Was ist Statistik?

EinbettenHerunterladen
Einfu¨hrung in die Statistik
mit Beispielen aus der Biologie
Thomas Fabbro
“The aim of computing is insight, not numbers.”
Was ist Statistik?
Die Statistik als Disziplin (“statistics”) besch¨aftigt sich mit dem
Sammeln, Organisieren, Analysieren, Interpretieren und
Pr¨asentieren von Daten (nach Dodge, Cox und Commenges 2006).
Wieso brauchen wir eigentlich Statistik?
Arten von Variablen
Messbare und Z¨ahlbare Variablen
numeric f¨
ur kontinuierlich Variablen, alle Zwischenschritte
sind m¨
oglich
integer f¨
ur Ganze Zahlen
Kategorielle Variablen
factor f¨
ur Kategorien (z. B. “Fabaceae”, “Rosaceae”,
“Apiaceae”).
logical Eine Variable die nur die Werte TRUE oder FALSE
annehmen kann (z. B. “m¨annlich”, “weiblich”).
Diese Einteilung basiert auf der Klasseneinteilung von R.
Beispiele in R
> weight <- c(0.001, 100, 3200000, 1000, 2.56, 0.001,
+
0.01)
> legs <- as.integer(c(0, 0, 4, 0, 0, 8, 6))
> kingdom <- factor(c("animal", "fungi", "animal",
+
"animal", "plant", "animal", "animal"))
> animal <- c(TRUE, FALSE, TRUE, TRUE, FALSE, TRUE,
+
TRUE)
Es gibt zwei Arten wie man Variablen beschreiben kann:
Kenngr¨
ossen
graphischen Darstellungen
Charakteristika kontinuierlicher Variablen
Lage
Streuung
Form
H¨aufung (”cluster”) Werte treten in Klumpen auf.
K¨
ornung (”granularity”) Nur bestimmte Werte treten auf.
Boxplot
Boxplot
●
0
1
2
3
4
gemessene Werte
5
6
Boxplot
Drawing the box:
Find the median. Then find the median of the data
values whose ranks are less than or equal to the rank of
the median. This will be a data value or it will be half
way between two data values.
Drawing the whiskers:
The maximum length of each whisker is 1.5 times the
interquartile range (IQR). To draw the whisker above the
3rd quartile, draw it to the largest data value that is less
than or equal to the value that is 1.5 IQRs above the 3rd
quartile. Any data value larger than that should be
marked as an outlier.
Der Boxplot wurde von Tukey eingef¨
uhrt. Heute gibt es viel
verschiedene Formen. Es ist daher gut, wenn man immer angibt
wie man den Boxplot konstruiert hat.
Histogramm
0.0 0.1 0.2 0.3 0.4
Dichte
Histogramm
0
1
2
3
4
gemessene Werte
5
6
Histogramm
Um ein Histogramm zeichnen zu k¨onnten muss man folgende zwei
Punkte festlegen:
K¨astchenbreite bzw. die Anzahl der K¨astchen (K¨astchen:
“bin”)
Startpunkt
Auch f¨
ur dieselben Werte sieht nicht jedes Histogramm gleich aus!
Histogramm
Anzahl Kästchen: 4
0.00
0.00
0.05
0.10
frequency
0.10
0.05
frequency
0.15
0.15
Anzahl Kästchen: 2
0
2
4
6
8
10
12
0
4
6
8
10
12
Anzahl Kästchen: 16
0.20
0.15
frequency
0.10
0.10
0.00
0.05
0.05
0.00
frequency
0.15
0.25
0.20
0.30
Anzahl Kästchen: 8
2
0
2
4
6
8
10
12
0
2
4
6
8
10
12
Histogramm
Startpunkt: −0.6
0.15
Density
0.00
0.00
0.05
0.10
0.10
0.05
Density
0.15
0.20
0.20
Startpunkt: 0
2
4
6
8
10
12
0
2
4
6
8
10
x
x
Startpunkt: −1.2
Startpunkt: −1.8
12
0.10
Density
0.05
0.10
0.00
0.00
0.05
Density
0.15
0.15
0.20
0
0
2
4
6
8
10
12
−2
0
2
4
6
8
10
Empirische Dichte
0.0 0.1 0.2 0.3 0.4
Dichte
"gaussian kernel"
0
1
2
3
4
gemessene Werte
5
6
12
Wahrscheinlichkeitsverteilung
µ
1
x¯ =
n
n
xi
i=1
Kenngr¨ossen fu¨r kontinuierliche Variablen
Lage Mittelwert, Median, Modus
Streuung Spannweite, Quartilsabstand, Varianz
Form Schiefe: z. B. rechtsschief = linkssteil, linksschief =
rechtssteil
W¨
olbung: steilgipflig, flachgipflig
weitere Begriffe: symmetrisch, unimodal, bimodal,
multimodal
Wahrscheinlichkeitsverteilung
Wahrscheinlichkeitsverteilung
0.25
χ2df=3
Dichte
0.20
6.25
0.15
0.10
0.05
0.1
0.04
0.00
0
1
2
3
4
5
x
6
7
8
9
Wahrscheinlichkeitsdichte
N (µ = 0, σ 2 = 1)
Dichte
0.4
1.96
0.3
0.2
0.1
0.025
−2
−1
0
1
x
Quantile-Quantile-Diagramm
Theoretische Quantile
(der Normalverteilung)
●
●
1
●
●
●
0
●
●
●
−1
●
●
0
1
2
3
4
5
Empirische Quantile
entstprechen den
geordneten Beobachtungen
6
2
2
2
2
●
1
●
●
●
0
●
●
●
−1
●
●
Theoretische Quantile
●
Theoretische Quantile
●
1
●
●
●
0
●
●
●
−1
●
●
2
4
6
0
2
4
0
●
●
●
0
●
●
●
●
1
●
●
●
●
0
●
●
●
−1
●
●
●
6
●
●
●
−1
●
●
0
2
4
6
0
Empirische Quantile
2
2
●
●
●
●
●
●
1
●
●
●
0
●
●
●
−1
●
●
●
●
2
4
6
●
0
●
●
●
−1
●
●
−2
−2
0
Empirische Quantile
●
1
●
−2
6
2
Theoretische Quantile
●
●
●
0
4
●
Theoretische Quantile
●
1
2
Empirische Quantile
●
0
●
0
−2
Empirische Quantile
−1
●
1
●
−2
4
6
2
●
●
−2
2
4
●
Theoretische Quantile
1
2
Empirische Quantile
2
●
0
●
6
●
Theoretische Quantile
●
●
−1
Empirische Quantile
2
−1
●
−2
Empirische Quantile
Theoretische Quantile
●
●
●
0
●
−2
0
●
1
●
−2
Theoretische Quantile
Theoretische Quantile
●
2
4
6
0
Empirische Quantile
2
4
6
Empirische Quantile
Verteilungsformen und Q-Q-Diagramme
density
normal quantiles
density
normal quantiles
skewed to the left
density
normal quantiles
skewed to the right
density
normal quantiles
platykurtic
leptokurtic
density
normal quantiles
bimodal
density
normal quantiles
normal
Verteilung
0.08
Dichte
0.06
0.04
0.02
0.00
150
160
170
180
Körpergrösse bei Frauen
Verteilung eines Mittelwertes
1.0
Dichte
0.8
0.6
0.4
0.2
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●●
●●
●●
●
●●
●
●
●
●●
●●
●
●
●
●
●
●
●
●
●
●
●●
●
●●
●
●
●
●
●●
●
●
●
●●
●
●
●●
●
●
●
●
●
●●
●
●●
●
●
●
●
●
●
●
●
●●
●
●●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
0.0
162
164
166
168
170
x von 300−Körpergrössen (Frauen)
172
Verteilung eines Mittelwertes
1.0
0.8
Dichte
0.6
0.4
0.2
0.0
162
164
166
168
170
172
x von 30 bzw. 300−Körpergrössen bei Frauen
Verteilung einer Statistik
s2
=
Varianz
n
¯)2
i=1 (xi − x
1
n−1
Dichte
Dichte
Mittelwert
x¯ = n1 ni=1 xi
x
s2
Jede Statistik folgt einer “eigenen” Verteilung
Die Streuung ist von der Gr¨osse, n, der Stichprobe abh¨angig
Die Form der Verteilung ist von der Gr¨osse der Stichprobe
unabh¨angig
Der Zentrale Grenzwertsatz
1.0
Dichte
0.8
0.6
0.4
0.2
2
4
6
8
10
Mittelwerte aus 50 Messwerten
Der Zentrale Grenzwertsatz
Die Verteilung des Mittelwertes aus n-Messwerten n¨ahert
sich f¨
ur wachsende n immer mehr einer Normalverteilung
und dies unabh¨angig von der Verteilung aus welcher die
Messwerte gezogen wurden.
Vertrauensintervall
1.0
Dichte
0.8
0.6
0.4
0.2
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
● ● ●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
● ●●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
0.0
162
164
166
168
170
172
x (n = 30) Mittlere Körpergrösse
Vertrauensintervall
Jetzt starten wir wie im richtigen Leben mit einer einzelnen
Stichprobe:
x1 = 166.1, x2 = 178.8, x3 = 169.5, x4 = 165.9, x5 = 172.1, x6 =
177.3, x7 = 165.2, x8 = 164.3, x9 = 175.6, x10 = 173.7, x11 =
171.8, x12 = 172.2, x13 = 172.6, x14 = 162.6, x15 = 168.6, x16 =
172, x17 = 161.3, x18 = 169.7, x19 = 160, x20 = 170.1, x21 =
165.1, x22 = 172.9, x23 = 168.1, x24 = 167.5, x25 = 180.1, x26 =
172.2, x27 = 157.8, x28 = 177.2, x29 = 167.4, x30 = 174.6
Vertrauensintervall
Zwei Wege:
Wir k¨
onnen eine Annahme Treffen u
¨ber die
Wahrscheinlichkeitsverteilung aus welcher wir die Stichprobe
gezogen haben. Dann k¨onnten wir beliebig oft eine Stichprobe
der selben Gr¨
osse ziehen, den Mittelwert berechnen und so die
Verteilung der Mittelwerte ermitteln.
Wir k¨
onnen den Zentralen Grenzwertsatz anwenden. Dazu
m¨
ussten wir aber die beiden Parameter der Normalverteilung
unseres Mittelwertes besser kennen, namentlich den
Mittelwert und die Varianz.
Vertrauensintervall: Die Lage
Mangels besserer Informationen w¨ahlen wir den Mittelwert unserer
Stichprobe (¯
x ) als Erwartungswert f¨
ur die Mittelwerte, x¯.
Vertrauensintervall: Die Streuung
Der Standardfehler ist die Standardabweichung einer Statistik.
Meistens spricht man vom Standardfehler und meint damit die
Standardabweichung des Mittelwertes, welche folgendermassen
berechnet wird:
sx
sx¯ = √
n
Vertrauensintervall
Mit diesen Angaben k¨
onnen wir nun aus unserer Stichprobe die
Verteilung des Mittelwertes veranschaulichen.
Dichte
0.3
0.2
0.1
0.0
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
● ● ●
●
●
●
●
●
●
●
●
●●
●●
●●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
● ●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●●
●
164
166
168
170
172
x (n = 30) Mittlere Körpergrösse
Die mittlere K¨
orpergr¨
osse betr¨agt 169.7, 95%-Vertrauensintervall:
[167.8; 171.7].
Vertrauensintervall
Die Formel f¨
ur das 1 − α-Vertrauensintervall f¨
ur den Mittelwert
lautet:
s
s
x¯ − z(1−α/2) √ ; x¯ + z(1−α/2) √ .
n
n
Mit zα bezeichnen wir das α-Quantil einer
Standardnormalverteilung, N (µ = 0, σ 2 = 1). Wichtig ist
besonders der Wert z(1−α/2) = 1.96.
Zum Merken:
Das 95%-Vertrauensintervall umfasst den Mittelwert
plus/minus zwei mal den Standardfehler.
Zusammenh¨ange
Gibt es einen Zusammenhang . . .
ein
Feuerzeug
mit sich
tragen
Lungenkrebs
Aspekte von Zusammenh¨angen
St¨arke
Folgerichtigkeit
Spezifit¨at
Zeitlichkeit
Biologischer Gradient
Plausibilit¨at
Stimmigkeit
Experiment
Analogie
Austin Bradford Hill, The Environment and Disease: Association or Causation? Proceedings of the Royal
Society of Medicine, 58 (1965), 295-300.
Zuf¨allige und systematische Fehler
Zuf¨alliger Fehler
Systematischer Fehler, Bias
Diese beiden Fehler treten in Kombination auf, zwei Beispiele:
Systematischer Fehler (Bias): Definition
“Any process at any stage of inference which tends to produce
results or conclusions that differ systematically from the truth.”
Nach Murphy, The Logic of Medicine, Baltimore: John Hopkins University Press, 1976 aus Sackett, Bias
in Analytic Research, Journal of Chronic Disease, 1979.
Wann kann ein systematischer Fehler auftreten?
Literatursuche
Festlegen und Ausw¨ahlen der Studienpopulation
Durchf¨
uhren der experimentellen Intervention (“Exposure”)
Messen von “Exposure” und “Outcome”
Analysieren der Daten
Interpretieren der Analyse
Publizieren der Resultate
Basierend auf dieser Einteilung hat Sackett einen Katalog von 35
systematischen Fehlern erstellt.
Sackett, Bias in Analytic Research, Journal of Chronic Disease, 1979.
Design
Zusammenh¨ange
Raucher
sein
ein
Feuerzeug
mit sich
tragen
Lungenkrebs
Bias-Einteilung nach Struktur
Confounding Bias Common causes of exposure and outcome.
(“Exposure” und “Outcome” haben eine
gemeinsame Ursache.)
Selection Bias Conditioning on commen effects of exposure
and outcome.
(Ein gemeinsamer Effekt von “Exposure” und
“Outcome” wird ber¨
ucksichtigt.)
Information Bias Systematisch fehlerhafte Information u
¨ber
“Exposure”, “Outcome” oder andere Variabeln
welche f¨
ur die Studie herangezogen werden.
z.B. Hern´
an and Robins, A Structural Approach to Selection Bias, Epidemiology, 2004.
Confounding Bias
gemeinsame
Ursache
"Outcome"
"Exposure"
Anpassen des Studientypen: Zuf¨allige “Exposure”-Zuweisung
erlaubt es, “Confounding” durch bekannte und unbekannte
Variablen zu verhindern.
Anpassen der Datenanalyse: Nur mit Fachwissen ist es
m¨
oglich, “Confounder” zu identifizieren, zu messen und mit
gewissen Annahmen in der Analyse zu ber¨
ucksichtigen.
Selection Bias
Tod vor
der Geburt
Herzmissbildung
Folsäure
Selection Bias
gemeinsamer Effekt
"berücksichtigt"
"Outcome"
"Exposure"
Kann auch in experimentellen Studien mit zuf¨alliger
“Exposure”-Zuweisung vorkommen (z.B. durch fehlende
Messwerte, “loss to follow-up”).
Information Bias
Systematisch fehlerhafte Information u
¨ber “Exposure”, “Outcome”
oder andere Variabeln welche f¨
ur die Studie herangezogen werden.
Definition ist unabh¨angig von der kausalen Struktur.
Wenn es sich um eine kategoriale Messgr¨osse handelt, spricht
man auch von “Misclassification Bias”.
Verschiedene Formen k¨onnen unterschiedlich eingeteilt werden
(z.B. “differential / non-differential”, “Exposure / Outcome /
Covariate Misclassification”).
Hypothesen
Wissenschaftler formulieren gest¨
utzt auf Beobachtungen und viel
Fachwissen Hypothesen.
Beispiel:
Hypothese: Alle Schw¨ane sind weiss.
Falsifikationismus
August Weismann, 1868 meinte, es
l¨asst sich eine wissenschaftliche Hypothese zwar niemals
erweisen, wohl aber, wenn sie falsch ist, widerlegen, und
es fragt sich deshalb, ob nicht Thatsachen beigebracht
werden k¨
onnen, welche mit einer der beiden Hypothesen
in unaufl¨
oslichem Widerspruch stehen und somit dieselbe
zu Fall bringen.
Hypothesen Test: Analogie zu anderen Testverfahren
Nullhypothese H0 : Kein Feuer
Kein Alarm: H0 akzeptieren
Alternativhypothese HA : Feuer
Alarm: H0 verwerfen, HA akzeptieren
H0 akzeptieren kein Alarm
H0 verwerfen Alarm
H0 wahr
kein Feuer
(1 − α)
Typ I Fehler, α
HA wahr
Feuer
Typ II Fehler, β
(1 − β)
1 − α: Vertrauensniveau, misst Vertrauen, dass kein Alarm auch
wirklich kein Feuer bedeutet (hohes Vertrauen: falscher Alarm
selten)
1 − β: Power: misst Wahrscheinlichkeit, dass Feuer auch Alarm
ausl¨
ost
Wu¨rfelbeispiel
Nullhypothese Die Maschine zeigt den Mittelwert von 4 W¨
urfen
mit einem 20-seitigen W¨
urfel.
Alternativhypothese Die Maschine zeigt nicht den Mittelwert von
4 W¨
urfen mit einem 20-seitigen W¨
urfel.
H0
H0
akzeptieren
verwerfen
H0 wahr
(1-α)
Typ I Fehler, α
HA wahr
Typ II Fehler, β
(1-β)
Konventionen des statistischen Testens:
Wir lehnen die Nullhypothese ab, wenn ein bestimmte
“extreme” Statistik beobachtet wird.
Wir bezeichnen etwas als “extrem”, wenn es mit einer
Wahrscheinlichkeit von < 5 % auftritt (α < 0.05).
Gibt es einen Zusammenhang
zwichen dem Geschlecht und der K¨orperl¨ange des Dreistachligen
Stichlings in der Bodenseeregion.
“Exposure” Geschlecht
“Outcome” K¨
orperl¨ange
Nullhypothese Es gibt keinen Unterschied in der K¨orperl¨ange
zwischen m¨annlichen und weiblichen Fischen.
H0 : µ♀ − µ♂ = 0
Alternativhypothese Es gibt einen Unterschied in der K¨orperl¨ange
zwischen m¨annlichen und weiblichen Fischen.
HA : µ♀ − µ♂ = δ
Teststatistik Differenz in der K¨orperl¨ange zwischen m¨annlichen
und weiblichen Fischen. x¯♀ − x¯♂
Hypothesentesten
H0 : µweibl. − µmännl. = 0
kritischer Wert
Dichte
kritischer Wert
H1 : µweibl. − µmännl. = δ
β
α 2
"P−Wert"/ 2
−10
α 2
"P−Wert"/ 2
−5
0
δ
5
10
15
xweibl. − xmännl.
Alternative Szenarien
Szenario 1:
Dichte
Gr¨ossere Stichprobe und daher weniger Streuung.
x0 − x1
20
Alternative Szenarien
Szenario 1 (Forts.):
Dichte
Bei gen¨
ugend grossen Stichproben k¨onnte man das
Signifikanzniveau senken, ohne wesentlich an Power zu verlieren.
x0 − x1
Alternative Szenarien
Szenario 2:
Dichte
Wenn der Unterschied zwischen den beiden Gruppen gen¨
ugend
gross ist, k¨
onnen relativ kleine Stichprobengr¨ossen ausreichen, um
einen signifikanten Unterschied zu zeigen.
x0 − x1
Signifikant, oder nicht?
Ein statistischer Test soll helfen zu entscheiden, ob die Daten mit
der Nullhypothese im Einklang stehen oder nicht. Wie viele andere
Tests liefert eine statistischer Test in erster Linie eine “ja/nein”
Antwort und dazu braucht es ein Signifikanzniveau, α, welches die
kritische Grenze darstellt.
Signifikant, oder nicht?
Der “P-Wert” ist ein verfeinertes Mass ob ein Test signifikant ist
oder nicht. Leider birgt der “P-Wert” einige Probleme:
Interpretation Der “P-Wert” wird sehr h¨aufig falsch interpretiert.
Auch die Bezeichnung “P”, welche meist f¨
ur
Wahrscheinlichkeiten verwendet wird, ist tr¨
ugerisch,
da sich die Wahrscheinlichkeit nur auf das
hypothetische Wiederholen des Experimentes bzw.
der Untersuchung bezieht.
Vereinheitlichung F¨
ur jeden statistischen Test kann man einen
“P-Wert” berechnen. Jeder statistische Test besteht
jedoch aus einer Nullhypothese und einer Teststatistik
und die Wahl dieses Paares beeinflusst den “P-Wert”.
Es besteht also die Gefahr, dass man dieser Zahl
mehr Beachtung schenkt, als sie es wert ist und sie
allzu sorgenlos ohne weitere Angaben nutzt.
Signifikant, oder nicht?
Ein Mittelweg zwischen “ja/nein” und dem “P-Wert”: Die
Sternchen-Konvention
Interpretation
Notation
P > 0.05
nicht signifikant
(n.s.)
0.05 ≥ P > 0.01
schwach signifikant
*
0.01 ≥ P > 0.001 stark signifikant
**
0.001 ≥ P
sehr stark signifikant ***
Relevant, oder nicht?
Wie gross muss ein Unterschied sein, dass er in einem
statistischen Test signifikant wird?
Was bedeutet es, wenn bei einem geplanten Experiment der
“P-Wert” sehr stark signifikant wird (P ≤ 0.001).
Was bedeutet es, wenn bei einem geplanten Experiment ein
statistischer Test nicht signifikant wird (P > 0.05).
Es ist unerl¨asslich das Resultat einer Untersuchung u
¨ber die
Sch¨atzer zu interpretieren! Die alleinige Aussage, dass ein
Testresultat signifikant ist erlaubt es nicht dieses zu interpretieren.
Genauso kann man nicht sagen, dass ein Unterschied der in einer
Untersuchung nicht signifikant wurde nicht relevant sei.
“The abscence of evidence is no evidence for abscence.”
Carl Sagan
Ein- oder zweiseitig Testen?
Um den kritischen Wert zu finden f¨
ur das Verwerfen der
Nullhypothese haben wir an beiden Enden der Verteilung jeweils
eine Fl¨ache von α/2 abgetrennt. Es stellt sich nat¨
urlich die Frage,
wieso wir nicht einfach α auf einer Seite abtrennen, da wir ja
sowieso erwarten, dass der Effekt in die bekannte Richtung geht.
Sollen wir also ein- oder zweiseitig testen?
Ganz einfach - zweiseitig! Weil wir ja nicht sicher sind, dass der
Effekt in diese Richtung geht. Einseitig Testen darf man nur, wenn
der Effekt aus z.B. physikalischen Gr¨
unden nur in eine Richtung
gehen kann.
Lage-Vergleich zwischen zwei Stichproben
Ein pragmatischer Ansatz:
Immer den Rangsummentest von Wilcoxon verwenden.
Nullhyothese H0 : Yg ,i F(i.i.d.)
Die Verteilung F kann eine beliebige Verteilung sein,
aber alle Messwerte Yg ,i m¨
ussen aus der selben
Verteilung sein und unabh¨angig voneinander.
Alternativhypothese HA : Yg =1,i F1 , Yg =2,i F2 wobei
F2 (x) = F1 (x − δ) und δ = 0
Teststatistik Die Teststatistik W wird wie folgt berechnet:
Rg ,i =
Rang (Yg ,i |Yg =1,i=1 . . . Yg =1,i=n1 , Yg =2,i=1 . . . Yg =2,i=n2 )
1
Ug =1 = i=n
i=1 Rg =1,i
der gr¨
ossere der beiden Werte
Wa = n1 n2 + n2 (n22 +1) − Ug =1 und Wb = n1 n2 − Ug =1
Gibt es einen Unterschied zwischen weiblichen und
m¨annlichen Stichlingen?
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
sex
f
f
f
f
m
m
f
f
m
m
m
m
f
m
m
bodySize
88.22
90.95
78.63
81.85
69.47
87.33
75.06
78.64
65.63
75.91
82.8
80.14
84.18
74.68
70.59
rank
29
30
16
22
3
28
10
17
1
13
23
21
26
9
5
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
sex
m
m
m
m
m
m
m
f
f
f
f
f
f
f
f
bodySize
75.55
78.04
73.82
70.57
69.38
74.55
79.62
75.89
79.3
74.14
83.92
80.07
77.64
86.51
83.28
Gibt es einen Unterschied zwischen weiblichen und
m¨annlichen Stichlingen?
Dichte
0.015
0.010
0.005
0.000
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
0
50
100
150
200
W
Wahrscheinlichkeitsdichte der Wilcoxon-Statistik
(n1 = n2 = 15) unter der Nullhypothese
rank
11
15
6
4
2
8
19
12
18
7
25
20
14
27
24
Gibt es einen Unterschied zwischen weiblichen und
m¨annlichen Stichlingen?
> wilcox.test(bodySize ~ sex, data = stickleback,
+
conf.int = TRUE)
Wilcoxon rank sum test
data: bodySize by sex
W = 177, p-value = 0.006568
alternative hypothesis: true location shift is not equal to 0
95 percent confidence interval:
1.73 9.92
sample estimates:
difference in location
5.88
Der Sch¨atzer “difference in location” entspricht dem Median
der Differenzen zwischen allen m¨oglichen Wertepaaren die
man zwischen den beiden Gruppen bilden kann.
Vergleich der Power zwischen Wilcoxon und T-Test
Wilcoxon Test
t Test
50
50
N=47
●
45
Total sample size
0
0.8
35
30
25
N=44
●
5
0.9
40
0.9
r=
we
5
Po
0.8
0.8
.95
0
0.9
r=
we
Po
.85
Total sample size
45
40
35
30
25
θ=1
θ=1
20
20
0.6
0.8
1.0
Delta, δ
1.2
1.4
0.6
0.8
1.0
1.2
1.4
Delta, δ
F¨
ur diese Abbildungen habe ich Messwerte aus einer
Normalverteilung N (µ = 0, σ = 1) bzw. N (µ = δ, σ = 1) gezogen
und mit den beiden Tests verglichen.
Multiples Testen
H0 : µweibl. − µmännl. = 0
kritischer Wert
Dichte
kritischer Wert
α 2
−10
α 2
−5
0
5
10
15
20
xweibl. − xmännl.
Multiples Testen
Anzahl
Tests
(k)
(k)
Wenn H0 wahr ist verwerfen wir H0
f¨alschlicherweise mit folgenden
Wahrscheinlichkeiten
Wahrscheinlichkeiten
α
nie
mindestens einmal
1
0.95
0.05
1
0.95
0.05
0.05
2
0.9025
0.0975
2
0.9025
0.0975
0.0253
3
0.8574
0.1426
4
0.8145
0.1855
k
k
(1 − α )
1 − (1 − α )k
3
0.8574
0.1426
0.017
4
0.8145
0.1855
0.0127
k
k
k
(1 − α )
1 − (1 − α )
ˇ
Die Dunn-Sid´ak Methode erlaubt das Signifikanzniveau der
einzelnen Tests (α ) so anzupassen, dass die Fehlerrate f¨
ur eine
ganze Familie von Testen auf dem Niveau von α gehalten werden
kann.
α = 1 − (1 − α )k
α = 1 − (1 − α)1/k
Lage-Vergleich zwischen mehrerer Stichproben
Auch um die Lage von mehr als zwei Stichproben zu testen gibt es
einen empfehlenswerten Test basierend auf den R¨angen der
Messwerte, den Kruskal-Wallis-Test. Er tested ob die Lage aller
Stichproben gleich ist (Nullhypothese).
> library(asuR)
> data(pea)
> kruskal.test(length ~ trt, data = pea)
Kruskal-Wallis rank sum test
data: length by trt
Kruskal-Wallis chi-squared = 38.4368, df = 4, p-value
= 9.105e-08
Wir m¨
ochten uns jedoch hier auch die parametrische Form etwas
genauer anschauen weil diese die Grundlage ist f¨
ur das Verst¨andnis
der meisten Formen des “Modellierens”. Die Methode welche es
erlaubt Lageparameter (im speziellen Fall Mittelwerte) aus
mehreren Stichproben zu vergleichen nennt man Varianzanalyse
oder kurz ANOVA (f¨
ur “analysis of variance”).
Formulieren eines linearen Modells
Einweg-ANOVA:
E[yij ]
yij
=
∼
µi
indep. N (µi , σ 2 )
mit 1, . . . , i Kategorien, jede mit 1, . . . , j Untersuchungseinheiten.
Das tolle an dieser Formulierung ist, dass sie uns erlaubt den
Erwartungswert der Kategorien und die Streuung der Messwerte
separat zu beschreiben.
Eine ANOVA Session mit R
>
>
>
>
+
>
+
>
>
>
>
library(asuR)
data(pea)
bwplot(~length | trt, layout = c(1, 5), data = pea)
trtSS <- sum(10 * (tapply(pea$length, pea$trt,
mean) - mean(pea$length))^2)
resSS <- sum((pea$length - rep(tapply(pea$length,
pea$trt, mean), each = 10))^2)
totSS <- sum((pea$length - mean(pea$length))^2)
model <- lm(length ~ trt, data = pea)
summary.aov(model)
summary(model)
Wahrscheinlichkeitsdichte der F Statistik
Dichte
0.6
0.4
0.2
0.0
0
2
4
6
8
10
F4, 45
Wahrscheinlichkeitsdichte der F Statistik
(df1 = 4, df2 = 45) unter der Nullhypothese
Formulieren eines linearen Modells
Einweg-ANOVA: am konkreten Beispiel
E[yij ]
=
µi
yij
=
∼
β0 + β1 x1 + β2 x2 + β3 x3 + β4 x4
indep. N (µi , σ 2 )
mit den Kategorien i = 1, 2, 3, 4, 5, jede mit den
Untersuchungseinheiten j = 1, . . . , 10.
Kontraste
Wenn man die Lage von mehr als zwei Stichproben beschreiben
m¨ochte, ist man meist nicht nur an der Nullhypothese die besagt
dass alle Lageparameter gleich sind (H0 : µ1 = µ2 = . . . = µk )
interessiert.
Oft m¨
ochte man auch einzelne Lageparameter oder
Linearkombinationen von Lageparametern vergleichen, man nennt
diese Vergleiche Kontraste, einige Beispiele:
H0 : µ1 − µ2 = 0 und H0 : µ1 − µ3 = 0 und H0 : µ2 − µ3 = 0
Alle paarweisen Vergleiche.
H0 : µ1 − 12 (µ2 + µ3 ) = 0, oder besser geschrieben
H0 : 1µ1 − 12 µ2 − 12 µ3 = 0 zum Testen ob der Mittelwert in der
Gruppe 1 sich vom Mittel der Gruppen 2 und 3 unterscheidet.
Kontraste
Wenn man die Nullhypothese (die Lage aller Stichproben ist
identisch) verwerfen konnte, dann kann man weitere Kontraste
untersuchen und testen. Allerdings darf man nicht beliebige und
auch nicht beliebig viele Kontraste untersuchen. Die Regeln sind:
Immer nur k − 1 Kontraste k¨onnen getestet werden, wenn
man k Stichproben untersucht.
Die Kontraste m¨
ussen orthogonal zueinander stehen, d.h. das
Produkt der Koeffizienten muss sich zu Null aufsummieren.
Ein Beispiel mit zwei Kontrasten:
1µ1 − 21 µ2 − 21 µ3 = 0
0µ1 − 1 µ2 + 1 µ3 = 0
(1 × 0) + (− 12 × −1) + (− 12 × +1) = 0
Eine ANOVA Session mit R
>
+
+
+
+
>
>
>
>
contr <- rbind(`control-sugar` = c(1, -1/4, -1/4,
-1/4, -1/4), `pure-mixed` = c(0, 1/3, 1/3,
-1, 1/3), `monosaccharides-disaccharides` = c(0,
1/2, 1/2, 0, -1), `gluc-fruc` = c(0, 1, -1,
0, 0))
ortho(contr)
contrasts(pea$trt) <- mancontr(contr)
model <- lm(length ~ trt, data = pea)
summary(model)
Kontraste: Tukey’s HSD
Student (WS Gosset) hat die Verteilung der t-Statistik entdeckt
f¨
ur zwei Stichproben welche sich nicht im Mittelwert
unterscheiden. Wenn es nun k-Stichproben gibt dann gibt es
k(k − 1)/2 paarweise Vergleiche mit einem zugeh¨origen t-Wert.
Tukey hat die Verteilung der gr¨ossten dieser t-Statistiken
beschrieben. Damit kann man alle paarweisen Vergleiche testen
und der Typ I Fehler wird nicht gr¨osser als α. Man nennt die
Methode auch Tukey’s HSD (Honest Significant Difference).
>
>
>
>
>
m0 <- aov(length ~ trt, data = pea)
(t0 <- TukeyHSD(m0))
par.old <- par(mar = c(4, 6, 2, 2))
plot(t0, las = 1)
par(par.old)
Lineare Modelle
Einweg ANOVA:
E[yij ]
=
=
yij ∼
mit i Kategorien,
µi
β0 + β1 x1 + . . . + βi−1 xi−1
indep. N (µ + αi−1 , σ 2 )
jede mit j Untersuchungseinheiten.
Regression:
E[yij ] = µ(xi ) = α + βxi = β0 + β1 xi
yij ∼ indep. N (α + βxi , σ 2 )
mit j gemessenen Untersuchungseinheiten an i Punkten auf der
kontinuierlichen x-Achse.
Regression
Es wurde der Gewichtsverlust von Tribolium confusum Larven nach
sechs Tagen ohne Futter bei verschiedenen Luftfeuchtigkeiten
gemessen.
12
●
Gewichtsverlust, mg
10
●
●
●
8
6
4
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
2
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
0
0
20
40
60
80
100
Relative Luftfeuchtigkeit
Regression
E[yij ]
yij
=
∼
µ(xi ) = α + βxi = β0 + β1 xi
indep. N (α + βxi , σ 2 )
> summary(lm(weight.loss ~ humidity, data = tribolium))
Call:
lm(formula = weight.loss ~ humidity, data = tribolium)
Residuals:
Min
1Q
-5.1830 -1.4155
Median
0.4941
3Q
1.3464
Max
5.3518
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 8.202155
0.436580 18.787 < 2e-16 ***
humidity
-0.044845
0.007421 -6.043 4.7e-08 ***
--Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 2.028 on 79 degrees of freedom
Multiple R-squared: 0.3162, Adjusted R-squared: 0.3075
F-statistic: 36.52 on 1 and 79 DF, p-value: 4.697e-08
Regression
Hypothesen Hat die erkl¨arende Variable (X ) einen linearen
Zusammenhang mit der Zielgr¨osse (Y ). Dazu testen
wir die Nullhypothese H0 : β = 0.
Vorhersage Ein Regressionsmodell erlaubt uns auch f¨
ur eine
gegebene erkl¨arende Variable (X ) eine Vorhersagen
f¨
ur die Zielgr¨osse (Y ) zu machen. Die
Wahrscheinlichkeitsverteilung f¨
ur die Vorhersage
2
ˆ
lautet N (ˆ
α + βxi , σ
ˆ ).
R 2 Sch¨atzt wie stark der lineare Zusammenhang
zwischen der erkl¨arenden Variable (X ) und der
Zielgr¨
osse (Y ) ist. Der Wert liegt immer zwischen 0
und 1. Leider ist R 2 sehr stark von der gew¨ahlten
Spannweite der erkl¨arenden Variablen abh¨angig und
kann somit nicht verwendet werden um Modelle von
verschiedenen Datens¨atzen zu vergleichen.
Voraussetzungen
¨
Das Uberpr¨
ufen der Voraussetzungen ist ein wesentlicher
Bestandteil einer Regressionsanalyse. Es geht einerseits darum die
Verletzungen aufzudecken und falsche R¨
uckschl¨
usse zu vermeiden
andererseits aber auch darum ein besseres Modell zu finden.
Residuen
Hebelarm
Die Streuung der Residuen ist konstant
Erwartungswert der Residuen ist null
Residuen sind normalverteilt
* Residuen sind unabh¨angig voneinander
Alle Messwerte haben denselben Einfluss
Streuung der Residuen
●
2
●
●
●
●
●
"Standartisierte" Residuen
●
1
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
−1
●
●
●
●
●
●
0
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
−2
●
●
●
●
4
5
6
7
8
"fitted values"
Der Erwartungswert der Residuen ist null
Die Streuung der Residuen ist konstant
Verteilung der Residuen
●
4
●
● ●
●●
●
●●
●
●
●●●
Residuen
2
●●●
●●●●●
●●●
●●
●●●●●
●
●●●
●●●●●
●
●
●
●●
●●
●
●●
●●●
●
●●
●●
●
●
●
●●
●●●●●
●
●●
●
0
−2
●●
●●
−4
●
●
●
●
−2
−1
0
Theoretische Quantile
Residuen sind normalverteilt
1
2
Ausreisser
Es gibt zwei Arten von Ausreissern
Y – Ausreisser Messwerte welche weit entfernt vom
Erwartungswert liegen (grosse Residuen) aber einen
kleinen Hebelarm haben.
X – Ausreisser Messwerte welche einen sehr grossen Hebelarm
(“leverage”) haben. Der Hebelarm gibt an was f¨
ur ein
Potenzial die einzelnen Werte haben das Modell zu
beeinflussen. Im Gegensatz zu Werten mit grossen
Residuen m¨
ussen Werte mit grossem Hebelarm nicht
zwingend problematisch sein.
Ausreisser: Potenzial-Residuen
Der Potenzial-Residuen-Plot zeigt welche Rolle einzelne Punkte im
“fitting” Prozess spielen.
●
●
●
●
Residuen
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
Potenzial
Alle Messwerte haben denselben Einfluss
Ausreisser: Potenzial-Residuen-Diagramm
Interpretation des Potenzial-Residuen-Diagramm (nur als grobe
Faustregel):
Gibt es Punkte mit grossem Potenzial und grossen Residuen
(Quadrant oben links) ist alles in Ordnung.
Gibt es Punkte mit grossem Potenzial und kleinen Residuen
lohnt es diese genauer zu untersuchen, sie k¨onnten (!) das
Modell stark verzerren.
Gibt es Punkte mit grossen Residuen aber kleinem Potenzial,
dann stimmt das Modell f¨
ur diese Punkte nicht. Das kann am
Modell liegen oder an den Punkten
Graphical Excellence
Graphical excellence is that which gives to the viewer the greatest number of ideas in the shortest time with
the least ink in the smallest space.
Edward R Tufte
Menschliche Wahrnehmung von graphischen Elementen
Rang
1
2
3
4
5
6
7
8
Aspekt
Position entlang einer gemeinsamen Achse
Position auf identischen aber nicht ausgerichteten Achsen
L¨ange
Winkel
Steigung
Fl¨ache
Volumen
Farbs¨attigung (schlecht sortierbar, gut diskriminierbar)
Cleveland 1985, Graphical Perception and Graphical Methods for Analyzing Scientific Data, Science, 229,
828-833.
Einige Ideen
Zeigen Sie die Messwerte (e.g. rug()), falls diese zu zahlreich
sind kann man auch eine zuf¨allige Auswahl zeigen.
Lassen Sie den Betrachter u
¨ber den Inhalt nachdenken statt
u
¨ber die Methoden.
Zeigen Sie was im Zentrum der Untersuchung steht. Zeigen
Sie z.B. die Differenz zwischen zwei Behandlungen mit einem
Vertrauensintervall statt die Mittelwerte der einzelnen
Gruppen.
Ob der Nullpunkt auf der Achse abgebildet sein soll oder nicht
muss sorgf¨altig u
¨berdacht werden.
Mittelwert-Differenz Abbildung
●
Männchen
●
Weibchen
70
75
80
85
10
15
●
Differenz
0
5
Körpergrösse von Stichlingen im Bodensee
Fehlerbalken
Fehlerbalken (“error bars”) k¨onnen Verschiedenes zeigen:
Standardabweichung der Daten
Standardabweichung einer Statistik (Standardfehler)
Vertrauensintervall
a
●●●●●●●●●●●●●●●●●●●●
b
●
●
●
●
●
data set
c
●
●
●
●
●
●●●●●●●●●●●●●●●●●●
d
● ● ● ●●
90
a
b
c
d
●●
● ● ● ● ●● ●● ●
●
●●
●
100
110
●
●
●
●
90
100
110
Alle vier Datensets haben gleich viele Beobachtungen, denselben
Mittelwert und dieselbe Standardabweichung.
Theoriefragen zu den bisherigen Inhalten
Was ist der Vorteil der gesch¨atzten empirischen Dichte
gegen¨
uber einem Histogramm?
Erkl¨aren Sie den Begriff Vertrauensintervall.
Was bedeutet der Begriff “Power” genau?
Was f¨
ur Angaben braucht es um das Resultat eines
statistischen Tests zu interpretieren?
Ein Beispiel: Sie haben getestet ob sich die K¨orpergr¨osse bei
Stichlingen hinsichtlich der Herkunft unterscheidet und dazu
je 100 Stichlinge aus zwei Zufl¨
ussen gemessen. Der
Unterschied in der Lage ist sehr stark signifikant. Welche
Angaben fehlen noch?
Sie m¨
ochten zeigen, dass die Stichlinge aus dem Obersee und
dem Untersee gleich gross sind. Wie gehen Sie vor?
Welchen Test verwenden Sie f¨
ur das Testen ob sich zwei
Stichproben hinsichtlich der Lage unterscheiden?
Theoriefragen zu den bisherigen Inhalten
Wie gross darf der P-Wert h¨ochstens sein, damit ein Resultat
relevant ist?
Was bedeutet beim statistischen Testen ein Typ I und II
Fehler?
Theoriefragen zu den bisherigen Inhalten
Was sind die drei wichtigsten Charakteristika kontinuierlicher
Variablen?
Wenn man die Wahrscheinlichkeitsdichte einer kontinuierlichen
Variablen aufzeichnet, welche Gr¨osse tr¨agt man auf der
y-Achse auf und wie interpretiert man diese Gr¨osse?
Beschreiben Sie was “Standardfehler” bedeutet.
Welche beiden Voraussetzungen m¨
ussen erf¨
ullt sein f¨
ur ein
Confounding?
Was ist der Unterschied zwischen einer Beobachtungsstudie
und einer experimentellen Studie?
Multivariate Statistik
Werden an einer Untersuchungseinheit mehrere Variablen
gemessen, dann erh¨alt man multivariate Daten.
Korrelation vs. Regression
Sowohl die Korrelation als auch die Regression untersuchen den
Zusammenhang zwischen zwei Variablen. Trotzdem unterscheiden
sie sich fundamental:
Korrelation untersucht den Zusammenhang zwischen zwei
“gleichberechtigten” Variablen (rY 1,Y 2 = rY 2,Y 1 ) und
ist ein rein beschreibendes Mass.
Regression untersucht den Einfluss von einer oder mehreren
“erkl¨arenden” Variablen, Xi , auf eine “Zielvariable”, Y .
Die Regressionsgerade von X auf Y ist in der Regel
eine andere als von Y auf X . Eine Regression erlaubt
nicht nur eine Beschreibung des Zusammenhangs,
sondern auch eine Vorhersage f¨
ur y − Werte wenn
man die entsprechenden x − Werte kennt.
Pearson Korrelation
Die Pearson Korrelation,
rY 1,Y 2
1
=
(n − 1)sy 1 sy 2
n
y 1i y 2i =
i=1
sy 1y 2
,
sy 1 sy 2
(1)
wobei sy 1 , sy 2 die Standardabweichung von y 1 und y 2 ist und sy 1y 2
die Kovarianz zwischen y 1 und y 2.
Die Pearson Korrelation liegt zwischen -1 und 1.
Die Pearson Korrelation misst “nur” den linearen Zusammenhang.
Achtung, auch ein sehr starker nicht linearer Zusammenhang kann
eine Pearson Korrelation rY 1,Y 2 = 0 haben!
Die Pearson Korrelation ist ganz und gar nicht robust und man
muss immer ein Auge darauf haben, ob es Werte gibt welche die
Korrelation stark beeinflussen (Scatterplot!).
Pearson Korrelation
> library(MASS)
> y1y2 <- mvrnorm(100, mu = c(0, 0), Sigma = matrix(c(1,
+
0.5, 0.5, 1), nrow = 2))
> cor(x = y1y2, method = "pearson")
[,1]
[,2]
[1,] 1.000000 0.531472
[2,] 0.531472 1.000000
Spearman Rangkorrelation
Die Spearman Rangkorrelation,
rY 1,Y 2 =
sRangY 1 ,RangY 2
,
sRangY 1 sRangY 2
(2)
misst im Gegensatz zur Pearson Korrelation nicht wie stark ein
“linearer”, sondern ein “monotoner” Zusammenhang zwischen den
Variablen Y 1 und Y 2 ist. Ein weiterer Vorteil der Spearman
Rangkorrelation ist die Robustheit.
> cor(x = y1y2, method = "spearman")
[,1]
[,2]
[1,] 1.0000000 0.5246205
[2,] 0.5246205 1.0000000
Interpretation einer Korrelation
Allgemein ist es schwierig Korrelationen zu interpretieren. Wichtig
zu beachten ist, dass ein gefundener Zusammenhang noch bei
weitem kein kausaler, urs¨achlicher Zusammenhang sein muss, auch
wenn der Zusammenhang sehr stark ist (siehe auch Design und
Bias).
Eine Korrelation ist ein Sch¨atzer f¨
ur einen Parameter - dazu
¨
braucht es eine Population. Die Uberlegung
wie diese Population
definiert werden kann ist oft hilfreich bei der Interpretation.
Wozu verwendet man Regressionen
Studium der Zusammenh¨ange (um die Variation der
Zielvariable erkl¨aren zu k¨onnen)
Vorhersage (um die Zielvariable m¨oglichst genau vorhersagen
zu k¨
onnen)
Statistische Kontrolle (um die Zielvariable m¨oglichst genau
einstellen zu k¨
onnen)
Multiple lineare Regression
Untersucht den Zusammenhang zwischen einer Zielgr¨osse und
mehreren erkl¨arenden Variablen.
E[yijk ]
yijk
=
∼
µ(xi ) = β0 + β1 x1i + β2 x2j
indep. N (β0 + β1 x1i + β2 x2j , σ 2 )
Tribolium Beispiel
ANCOVA
Untersucht den Zusammenhang zwischen einer Zielgr¨osse und
mindestens einer kontinuierlichen und einer kategoriellen
erkl¨arenden Variable.
E[yijk ]
yijk
=
∼
β0 + β1j + β2 xi
indep. N (β0 + β1j + β2 xi , σ 2 )
ANCOVA mit Interaktion
Wenn der Einfluss einer erkl¨arenden Variable von einer weiteren
Variable abh¨angt, spricht von von einer Interaktion dieser beiden
Variablen.
F¨
ur eine ANCOVA w¨
urde das folgendermassen aussehen:
E[yijk ]
yijk
=
∼
β0 + β1j + β2 xi + β3 xi j
indep. N (β0 + β1j + β2 xi + β3 xi j, σ 2 )
¨
Eine Ubung:
Fitten Sie eine Regression um mit den Variablen “Gr¨osse” und
“Geschlecht” das “Gewicht” zu erkl¨aren.
Gr¨
osse
150
165
180
170
Geschlecht
weiblich
m¨annlich
m¨annlich
weiblich
Gewicht
50
55
65
70
Zweiweg-ANOVA
E[yijk ]
yijk
=
∼
β0 + β1i + β2j
indep. N (β0 + β1i + β2j , σ 2 )
mit 1, . . . , i Kategorien im ersten Faktor und 1, . . . , j Kategorien
im zweiten Faktor, jeweils mit 1, . . . , k Untersuchungseinheit.
Zweiweg-ANOVA mit Interaktion
E[yijk ]
yijk
=
∼
β0 + β1i + β2j + β3ij
indep. N (β0 + β1i + β2j + β3ij , σ 2 )
Mit diesem Modell kann der Einfluss vom ersten Faktor auf den
Outcome abh¨angig vom zweiten Faktor modelliert werden.
90
Interaktions-Plot
60
70
G
g
40
50
Ertrag
80
Allel G
a
A
Allel A
Datentransformation
Warum transformieren wir?
Damit wir die gesch¨atzten Werte besser interpretieren k¨onnen
(Einheiten, Zentrieren, Skalieren).
Damit ein Zusammenhang linear wird
(und wir ihn mit den bekannten Methoden modellieren
k¨
onnen).
Damit die Residuen “normalverteilt” und die Variation
konstant wird
(und daher die Annahmen unserer Methoden nicht verletzt
werden).
Exkurs: Allometrie
Beziehung zwischen der Gr¨osse und der Form, Anatomie,
Physiologie und auch des Verhaltens eines Organismus.
Die klassische Formel
y = ax b ,
wird in der logarithmischen Form zu
log y = log a + b log x.
In dieser Form k¨
onnen wir u.a. durch Anpassen eines linearen
Modells die Parameter sch¨atzen. Wenn eine Skalierung isometrisch
ist erhalten wir eine Steigung von b = 1. Was f¨
ur eine Steigung
erwarten wir bei einer isometrischen Skalierung, wenn y eine
L¨angenangabe ist und x ein Volumen (oder Gewicht) darstellt?
Durch Umformen ergibt sich log l = ... + b log l 3 = ... + 3b log l
und daher eine Steigung von 3.
Linearit¨at erreichen
y = ax b ,
wird zu log y = log a + b log x.
4
4
b>1
3
3
y
y
b=1
2
0<0<1
1
2
b < −1
1
0
0
1
2
x
3
4
b = −1
−1 < b < 0
0
0
1
2
x
3
4
Linearit¨at erreichen
y = ae bx ,
wird durch logarithmieren log y = log a + bx.
1.0
b>0
50
b<0
0.8
0.6
30
y
y
40
0.4
20
10
0.2
0
0.0
0
1
2
3
4
0
1
x
Datentransformation: Logarithmieren
0.4
0.2
0.0
Dichte
0.6
"untransformiert"
0
2
4
6
8
10
gemessene Werte
0.0
0.2
0.4
0.6
"log−transformiert"
Dichte
2
x
−2
−1
0
1
log(gemessene Werte)
2
3
4
Datentransformation: Faustregeln
log(x)
log(x + 1)
log(10000x + 1)
x >1
x ≥0
1>x ≥0
√
f¨
ur Z¨ahldaten (Poisson)
x
x+
1
2
falls mit Null
√
arcsin x
1
1+x
2 ln 1−x
Potenztransformationen:
xn
√
x
ln(x)
√
1/ x
1/x
1/x n
0 < x < 1, Prozente und Proportionen
−1 < x < 1, Regressionskoeffizienten
linksschief
schwach rechtsschief
...
...
...
starkt rechtsschief
Mosteller und Tukey’s bulging rule (W¨olbungsregel)
Y
X
X
Y
Voraussetzungen
E[yij ]
yij
=
∼
µ(xi ) = α + βxi = β0 + β1 xi
indep. N (α + βxi , σ 2 )
Die Streuung der Residuen ist konstant
Erwartungswert der Residuen ist null
Residuen sind normalverteilt
Alle Messwerte haben denselben Einfluss/Hebelarm.
Residuen sind unabh¨angig voneinander
Körpergewicht (in kg)
Nicht unabh¨angige Residuen
120
120
●
110
110
●
100
90
90 ●
80
80 ●
GruppeA
Nicht unabh¨angige Residuen
> summary(lm(y ~ 1))$coefficients
(Intercept)
Estimate Std. Error t value
Pr(>|t|)
100
9.128709 10.95445 0.001628625
Körpergewicht (in kg)
Nicht unabh¨angige Residuen
120
120
●
110
110
●
●
122
●
107
100
90
90 ●
80
80 ●
●
92
●
79
GruppeA
Nicht unabh¨angige Residuen
> summary(lm(y2 ~ 1))$coefficients
(Intercept)
Estimate Std. Error t value
Pr(>|t|)
100
6.032649 16.57647 7.101236e-07
Nicht unabh¨angige Residuen
> library(geepack)
> summary(geeglm(y2 ~ 1, id = id, data = repData))$coefficients
(Intercept)
Estimate Std.err
Wald Pr(>|W|)
100 7.962804 157.7132
0
Was fu¨hrt zu nicht unabh¨angigen Residuen?
Information in den Daten welche nicht im Modell ber¨
ucksichtigt
wird.
Zeitliche Abh¨angigkeiten, z.B. Saisonalit¨at
R¨aumliche Abh¨angigkeiten
Genetische Abh¨angigkeiten
...
Was sind die Folgen von nicht unabh¨angigen Residuen?
Die Varianz welche gesch¨atzt wird und sich u.a. auf die
Standardfehler der Sch¨atzer und die Signifikanztests u
¨bertr¨agt ist
falsch.
Was gibt es fu¨r L¨osungen
Es muss eine andere Form der Modellierung verwendet werden
welche es erlaubt die Struktur der Daten richtig zu ber¨
ucksichtigen.
Mixed Effects Modelle
Generalized Estimating Equations (GEE)
...
Ein Beispiel
Sie nehmen auf zwei Feldern (“Field”) jeweils drei Bodenproben
(beschriftet mit “MeasurementID”) und bestimmen den
Stickstoffgehalt. Nun m¨
ochten Sie diese Daten analysieren.
> library(geepack)
> m.gee <- geeglm(Soilnitrogen ~ Field, id = MeasurementID,
+
data = na.omit(nitrogen))
Analyse von kategoriellen Daten
Was f¨
ur Arten von Variabeln kennen Sie?
Welche Arten von kategoriellen Daten kennen Sie?
Die wichtigsten Verteilungen von kategoriellen Daten
Binomial Verteilung
Multinominal Verteilung
Poisson Verteilung
Binomial Verteilung
Wir gehen davon aus, dass es n unabh¨angige und
identische Ziehungen y1 , y2 , ..., yn gibt.
Jede Ziehung ist entweder yi = 0 oder 1 (Misserfolg
/ Erfolg, FALSE / TRUE, “kein Event” / “Event”).
Die Wahrscheinlichkeit, dass eine Ziehung ein Erfolg
ist, bezeichnen wir als π.
Die Summer der Ziehungen folgt einer Binomial
Verteilung.
Mittelwert nπ
Varianz nπ(1 − π)
Identisch bedeutet, dass die Wahrscheinlichkeit eines Erfolges
bei allen Ziehungen dieselbe ist.
Unabh¨angig bedeutet, dass die Ziehung nicht von anderen
Ziehungen abh¨angt.
Multinominale Verteilung
Wie die binomial Verteilung aber jede unabh¨angige und identische
Ziehung f¨allt in eine von c Kategorien.
Poisson Verteilung
●
●
●
●
●
0.10
Häufigkeit
0.20
Wenn die Anzahl der Ziehungen nicht festgelegt ist wie bei der
binominal Verteilung, dann folgt die Anzahl der Erfolge einer
Poisson Verteilung. Beispiel: Anzahl der t¨otlichen Verkehrsunf¨alle
in Italien an einem Tag.
Eine Besonderheit der Poisson Verteilung ist, dass die Varianz dem
Mittelwert entspricht. Wenn der Mittelwert gross ist, dann n¨ahert
sich die Verteilung einer Normalverteilung.
Mittelwert = Varianz = µ
●
●
●
●
●
●
0.00
●
●
●
●
●
●
●
●
●
●
●
●
●
5
●
●
10
●
●
●
●
●
15
●
●
●
●
●
●
●
●
20
Overdispersal
H¨aufig beobachtet man Z¨ahldaten deren Streuung gr¨osser ist als
man aufgrund der Binomial oder Poisson Verteilung erwarten
w¨
urde. Dieses “Ph¨anomen” nennt man Overdispersal. Es kommt
daher, dass die einzelnen Ziehungen nicht “identisch” sind, d.h. aus
unterschiedlichen Population mit unterschiedlichen
Erwartungswerten kommen. Streng genommen handelt es sich
dann gar nicht um eine Zufallsvariable, sondern um eine Mischung
aus Zufallsvariablen.
Kreuztabellen, “contingency tables”
Familie
Fabaceae
Rosaceae
Wuchsform
Kraut Busch Baum
218
76
294
56
76
132
6
28
34
280
180
460
Beides (Familie und Wuchsform) sind Zielvariablen
(“response”, “outcome”) und nicht erkl¨arende Variablen.
Gibt es einen Zusammenhang zwischen Wuchsform und
Familie?
Der χ2 -Test oder der Fisher-Exakt-Test (nur f¨
ur 2 x 2
Tabellen) testen die
Nullhypothese: Es gibt keinen Zusammenhang.
Viel wichtiger als der Test ist allerdings die St¨arke des
Zusammenhangs.
Testen des Zusammenhangs
Pearson's Chi-squared test
data: table(plants$family, plants$type)
X-squared = 67.2916, df = 2, p-value = 2.442e-15
Was sagt uns dieser Test?
Was sagt uns dieser Test nicht?
Testen des Zusammenhangs
Familie
Kraut
Fabaceae
Rosaceae
218
76
294
Wuchsform
Busch
Baum
56
76
132
6
28
34
280
180
460
Wir erwarten Fabaceae und Rosaceae mit einer Proportion
von 0.61 und 0.39.
Wir erwarten die Wuchsform Kraut, Busch oder Baum mit
einer Proportion von 0.64, 0.29, 0.07.
Daraus k¨
onnen wir die Proportion f¨
ur die einzelnen Zellen
berechnen (“Erwartungswerte”).
Familie
Kraut
Wuchsform
Busch
Fabaceae
0.39
0.17
0.04
Rosaceae
0.25
0.11
0.03
Baum
Beschreibung von Zusammenh¨angen kategorieller Variablen
Daten aus einem Pestizidexperimet mit Helothis virescens, einem
Wuchsform leben.
Eulenfalter dessen Larven auf Tabakpflanzen
Familie
Kraut
Busch
Baum
¨
Uberleben
Geschlecht
tot 0.17lebend
Fabaceae
0.39
0.04
(0.47)
(0.12)
(0.01)
m¨annlich
4
16
20
weiblich
2
200.03
Rosaceae
0.25
0.1118
34
40(0.06)
(0.17) 6
(0.17)
Differenz zwischen zwei Proportionen, π1 − π2
Wuchsform
Relative
Risk,
π
/π
1
2
Familie
Kraut
Busch
Baum
Odds Ratio
Fabaceae 0.39
0.17
0.04
(0.47)
(0.12)
(0.01)
[2.92]
[-2.72]
[-3.23]
Rosaceae 0.25
0.11
0.03
(0.17)
(0.17)
(0.06)
[-3.64]
[3.39]
[4.03]
Um den Ursprung des Testresultates zu untersuchen
berechnen wir die “Erwartungswerte” f¨
ur die einzelnen Zellen
und vergleichen diese mit den Beobachtungen; “Nomenklatur”
Differenz zwischen zwei Proportionen, π1 − π2
Wird auch als “absolute risk difference” bezeichnet.
Wickler Beispiel: pmale − pfemale = 4/20 − 2/20 = 1/10
> prop.test(x = c(4, 2), n = c(20, 20))
2-sample test for equality of proportions with
continuity correction
data: c(4, 2) out of c(20, 20)
X-squared = 0.1961, df = 1, p-value = 0.6579
alternative hypothesis: two.sided
95 percent confidence interval:
-0.1691306 0.3691306
sample estimates:
prop 1 prop 2
0.2
0.1
Wertebereich: −1, . . . , 1
Relative Risk, π1 /π2
Es ist m¨
oglich, dass der Unterschied zwischen zwei Differenzen
anders zu interpretieren ist, wenn beide π klein oder gross sind, als
wenn beide nahe bei 0.5 sind.
Wickler Beispiel: pmale /pfemale = (4/20)/(2/20) = 2
Das 95% Vertrauensintervall dazu liegt bei [0.43; 7.47].
Wertebereich: 0, . . . , ∞
Odds Ratio
F¨
ur ein Ereignis mit der Wahrscheinlichkeit π sind die Odds
definiert als π/(1 − π), wie beim Wetten. Im Beispiel sind also die
Odds f¨
ur die M¨annchen (4/20)/(16/20) = 0.25. Das bedeutet,
dass es f¨
ur ein M¨annchen vier mal weniger wahrscheinlich ist zu
sterben, als zu u
¨berleben. Das Verh¨altnis
θ=
π1 /(1 − π1 )
π2 /(1 − π2 )
nennt man “odds ratio”.
Es betr¨agt in unserem Beispiel
((4/20)/(16/20))/((2/20)/(18/20)) = 2.25
Das Odds Ratio ¨andert nicht, wenn wir Zeilen und Spalten
vertauschen, d.h. wir m¨
ussen uns nicht auf eine Zielvariable
festlegen.
Wertebereich: 0, . . . , ∞
Struktur von glm’s (generalized linear model)
E[yij ]
yij
=
∼
β0 + β1 x1 + β2 x2 + . . . + βn xn
indep. N (β0 + β1 x1 + β2 x2 + . . . + βn xn , σ 2 )
Diese Form kennen wir von Regressionen und in fast
identischer Form kennen wir die Form auch von ANOVAs und
ANCOVAs.
Das Anwendungsspektrum l¨asst sich jedoch stark erweitern,
wenn man zus¨atzliche Flexibilit¨at einf¨
uhrt.
- andere Verteilungen
- eine Funktion welche den Zusammenhang zwischen dem
“linearen predictor” und dem Erwartungswert E[yij ] beschreibt.
E[yi ] = h(ηi )
Die Funktion h bezeichnen wir als “response function”. H¨aufig
wird jedoch die Umkehrfunktion h−1 = g beschrieben, welche
als “link function” oder einfach als “link” bezeichnet wird.
Struktur von glm’s (generalized linear model)
F¨
ur eine “gew¨
ohnliche” Regression w¨
urde also
lauten:
ηi = β0 + β1 x1 + β2 x2 + . . . + βn xn
E[yi ] = h(ηi ) = ηi
yi ∼ indep. N (ηi , σ 2 )
die Beschreibung
“linear predictor”
“response function”
“Verteilungsannahme”
Logistische Regression
Ein glm welches sehr h¨aufig verwendet wird ist die logistische
Regression.
ηi = β0 + β1 x1 + β2 x2 + . . . + βn xn
eη
E[yi ] = πi = h(ηi ) = 1+e
η
yi ∼ indep. B(πi )
F¨
ur den “linar predictor” ηi gilt das selbe wie bei einer
Regression bzw. ANOVA.
Als “link function” w¨ahlen wir die sogennante “logit” Funktion,
x
) (Der Name “logit” stammt von log
h−1 (x) = g (x) = log( 1−x
wie Logarithmus und Odds).
Als Verteilung w¨ahlt man eine Bernoulli Verteilung, B(πi ), mit
einem Parameter, der Wahrscheinlichkeit π des Ereignisses.
Logistische Regression, Interpretation der Koeffizienten
Wir betrachten einige Umformungen:
πi
= α + βxi
log 1−π
i
πi
= e α+βxi
1−πi
= e α e βxi
Wenn wir nun x um eins erh¨ohen, d.h. x ersetzten durch x + 1,
πi
= e α e β(xi +1)
1−πi
= e α e βxi e β
Daraus sehen wir, dass sich die Odds um den Faktor e β ver¨andern,
wenn wir von x nach x + 1 gehen.
Logistische Regression mit kategoriellen Pr¨adiktoren
Geschlecht
tot
m¨annlich
weiblich
4
2
6
¨
Uberleben
lebend
16
18
34
20
20
40
Odds f¨
urs Sterben der M¨annchen: om =
4/20
(1−4/20)
Odds f¨
urs Sterben der Weibchen: of =
Odds Ratio: OR =
om
of
= 2.25
Daraus ergibt sich z.B. of OR = om
2/20
(1−2/20)
= 0.25.
= 0.111.
>
>
>
>
y <- cbind(dead = c(2, 4), alive = c(18, 16))
sex <- factor(c("female", "male"))
glm1 <- glm(y ~ sex, family = "binomial")
summary(glm1)$coefficients
Estimate Std. Error
z value
Pr(>|z|)
(Intercept) -2.1972246
0.745356 -2.9478861 0.003199549
sexmale
0.8109302
0.931695 0.8703816 0.384091875
e αˆ = e −2.2 = 0.11 sch¨atzt uns die Odds f¨
ur die Weibchen mit
einem Standardfehler.
ˆ
e β = e 0.81 = 2.25 sch¨atzt uns das Oddsratio, d.h. den Faktor
mit welchem wir die Odds der Weibchen multiplizieren m¨
ussen
um die Odds der M¨annchen zu erhalten. Auch das Oddsratio
wird mit einem Standardfehler gesch¨atzt.
Logistische Regression mit kontinuierlichen Pr¨adiktoren
Im selben Versuch mit Heliothis virescens wurden 6 verschiedene
Dosen des selben Wirkstoffes untersucht.
Dosis
1 2 4
Geschlecht
m¨annlich
weiblich
1
0
4
2
9
6
8
16
32
13
10
18
12
20
16
Wir fokusieren uns zuerst mal nur auf die M¨annchen.
Logistische Regression mit kontinuierlichen Pr¨adiktoren
Eine m¨
ogliche Forschungsfrage k¨onnte lauten:
Nimmt die Sterblichkeit mit zunehmender Dosis signifikant zu?
Da wir erwarten, dass sich das Oddsratio jeweils gleich ver¨andert,
wenn wir die Dosis verdoppeln, werden wir den Logarithmus der
Dosis (mit der Basis 2) als Pr¨adiktor verwenden.
>
>
>
+
>
library(asuR)
data(budworm)
glm2 <- glm(cbind(num.dead, num.alive) ~ log2(dose),
data = budworm[budworm$sex == "male", ], family = "binomial")
summary(glm2)$coefficients
Estimate Std. Error
z value
Pr(>|z|)
(Intercept) -2.818555 0.5479868 -5.143472 2.697066e-07
log2(dose)
1.258949 0.2120655 5.936607 2.909816e-09
> summary(glm2)$coefficients
Wahrscheinlichkeit zu Sterben
Estimate Std. Error
z value
Pr(>|z|)
(Intercept) -2.818555 0.5479868 -5.143472 2.697066e-07
log2(dose)
1.258949 0.2120655 5.936607 2.909816e-09
1 4
9
13
18
1
20
●
●
Inf
9
0.8
●
0.6
1.86
●0.82
●
0.4
0.2
●
●
0
0.25
0.05
19 16
11
7
12 4
8
2
16
Dosis
0
32
e β = e 1.26 = 3.52 gibt an um welchen Faktor die Odds
zunehmen wenn wir die Dosis um den Faktor 2 erh¨ohen.
−α/β gibt uns die Dosis an wo die H¨alfte der Tiere sterben.
Multivariate logistische Regression
Selbstverst¨andlich w¨
urde es uns jetzt interessieren mit dem
Datensatz ein multivariates Modell zu rechnen. Damit k¨onnten wir
die Sterblichkeit mit zwei Pr¨adiktoren Dosis und Geschlecht
gleichzeitig erkl¨aren. Dadurch l¨asst sich u.a. die Frage beantworten
ob die Sterblichkeit unterschiedlich ist zwischen den Geschlechtern.
Multivariate logistische Regression
> glm3 <- glm(cbind(num.dead, num.alive) ~ log2(dose) +
+
sex, data = budworm, family = "binomial")
> summary(glm3)$coefficients
Wahrscheinlichkeit zu Sterben
Estimate Std. Error
z value
Pr(>|z|)
(Intercept) -3.473155 0.4685202 -7.413032 1.234445e-13
log2(dose)
1.064214 0.1310775 8.118971 4.701542e-16
sexmale
1.100743 0.3558271 3.093478 1.978249e-03
1
1 4
0 2
9
6
13
10
18
12
20
16
●
●
0.8
●
●
0.6
●
●
●
0.4
●
0.2
female
male
●
●
●
0
●
19 16
20 18
11
14
7
10
12 4
8
2
8
16
Dosis
0
4
32
... mit Interaktion
Ganz analog zu einer ANCOVA k¨onnen wir auch hier untersuchen,
ob die Zusahme der Sterbilichkeit mit steigender Dosis bei beiden
Geschlechtern gleich stark ist. Dazu fitten wir ein Modell mit einem
Interaktionsterm und testen ob dieser statistisch signifikant ist.
Multivariate logistische Regression
> glm4 <- glm(cbind(num.dead, num.alive) ~ log2(dose) *
+
sex, data = budworm, family = "binomial")
> round(summary(glm4)$coefficients, 3)
(Intercept)
log2(dose)
sexmale
log2(dose):sexmale
Estimate Std. Error z value Pr(>|z|)
-2.994
0.553 -5.416
0.000
0.906
0.167
5.422
0.000
0.175
0.778
0.225
0.822
0.353
0.270
1.307
0.191
Multivariate logistische Regression
> glm5 <- glm(cbind(num.dead, num.alive) ~ I(log2(dose) +
3) * sex, data = budworm, family = "binomial")
> round(summary(glm5)$coefficients, 3)
Estimate Std. Error z value
-0.275
0.231 -1.195
0.906
0.167
5.422
1.234
0.377
3.273
0.353
0.270
1.307
Pr(>|z|)
(Intercept)
0.232
I(log2(dose) - 3)
0.000
sexmale
0.001
I(log2(dose) - 3):sexmale
0.191
(Intercept)
I(log2(dose) - 3)
sexmale
I(log2(dose) - 3):sexmale
Multivariate logistische Regression
> glm6 <- glm(cbind(num.dead, num.alive) ~ I(log2(dose) +
3) + sex, data = budworm, family = "binomial")
> round(summary(glm6)$coefficients, 3)
(Intercept)
I(log2(dose) - 3)
sexmale
Estimate Std. Error z value Pr(>|z|)
-0.281
0.243 -1.154
0.249
1.064
0.131
8.119
0.000
1.101
0.356
3.093
0.002
Overdispersion
Wenn wir eine logistische Regression rechnen ist der sog.
“dispersion parameter” immer auf eins gesetzt. Das ist fast ein
bisschen so wie wenn wir in einer gew¨ohnlichen Regression den
Standardfehler, σ, auf eins setzten w¨
urden.
Um zu u
ufen, ob die beobachtete Variation gr¨osser ist, als wir
¨berpr¨
unter idealen Bedingungen einer logistischen Regression annahmen,
k¨
onnen wir den “dispersion parameter” sch¨atzten.
Dazu fitten wir ein quasibinomiales GLM.
Wir sehen wie der “dispersion parameter” gesch¨atzt wird.
Wir sehen ob das Ber¨
ucksichtigen des “dispersion parameters”
unsere Resultate qualitativ ver¨andert.
Overdispersion
> glm6.quasi <- glm(cbind(num.dead, num.alive) ~
+
I(log2(dose) - 3) + sex, data = budworm, family = "quasibinomial")
> summary(glm6.quasi)
Call:
glm(formula = cbind(num.dead, num.alive) ~ I(log2(dose) - 3) +
sex, family = "quasibinomial", data = budworm)
Deviance Residuals:
Min
1Q
Median
-1.10540 -0.65343 -0.02225
3Q
0.48471
Max
1.42944
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept)
-0.2805
0.1867 -1.503 0.16714
I(log2(dose) - 3)
1.0642
0.1006 10.574 2.24e-06 ***
sexmale
1.1007
0.2732
4.029 0.00298 **
--Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for quasibinomial family taken to be 0.5895578)
Null deviance: 124.8756
Residual deviance:
6.7571
AIC: NA
on 11
on 9
degrees of freedom
degrees of freedom
Number of Fisher Scoring iterations: 4
GLM Familien
Verteilung
gaussian
binomial
Gamma
poisson
inverse.gaussian
m¨
ogliche Link-Funktion
identity
log
inverse
logit
probit
cauchit
log
cloglog
inverse
identity
log
log
identity
sqrt
1/mu^2
inverse
identity
log
(logistisch)
(normal)
complementary log-log
Hinweis: Modelle mit dem selben linearen Pr¨adiktoren, η, aber mit
einer unterschiedlichen Link-Funktion, g , k¨onnen nur informell
verglichen aber nicht “getestet” werden.
Mixed Effects Models
Eine Unterscheidung welche nur auf kategorielle Pr¨adiktoren
zutrifft.
“Fixed Effects” Wenn es eine endliche Zahl an Kategorien (“levels”)
gibt denen wir einen Effekt zuordnen und wir an den
Effekten der einzelnen Kategorien interessiert sind.
Das Interesse liegt beim Sch¨atzen der Mittelwerte.
“Random Effects” Wenn es unendlich viele Kategorien gibt aber
nur ein zuf¨alliger Teil davon in unserer Stichprobe
gelandet ist und wir nicht an den Effekten der
einzelnen Kategorien interessiert sind. Das Interesse
liegt beim Sch¨atzen der Varianz.
Oft ist es direkt aus der Fragestellung m¨oglich zu unterscheiden ob
ein “Random Effect” oder ein “Fixed Effect” vorliegt.
Mixed Effects Models: Ein Beispiel
Sie m¨
ochten wissen, ob Knaben (5. Schuljahr) bessere Schulnoten
im Fach Mathe haben als M¨adchen.
Dazu haben Sie die Schulnoten eines ganzen Bezirkes erhalten und
wissen zu jedem Sch¨
uler das Geschlecht, die Bezeichnung der
Schule und die Bezeichnung der Klasse.
Schulnote, y Response
Geschlecht, β Fixed Effect
Schule, s Random Effect
Klasse, c Random Effect (genested in Schule)
E[ytijk | si , cij ]
ytijk | si , cij
=
∼
µ + βt + si + cij
indep. N (µ + βt + si + cij ;
σs2 + σc2 + σ 2 )
Bootstrapping
µ
1
x¯ =
n
n
xi
i=1
Bootstrapping
Name: “bootstrap” (Stiefelschlaufe), nach dem englischen
Spruch “pull oneself over a fence by one’s bootstraps”.
Eine Resamplingmethode zum bestimmen der Genauigkeit
eines Sch¨atzers.
Eine sehr einfache Methode die oft auch gut funktioniert.
Leider gibt es keine Garantie daf¨
ur!
Bootstrapping
Es gibt viele verschiedene Arten wie wiederholt eine Stichprobe
gezogen werden kann. Ich m¨ochte kurz auf das parametrische
und das nicht-parametrische Bootstrapping eingehen.
Das wiederholte Ziehen (“resampling”) f¨
uhrt immer dazu, dass
man neben dem Sch¨atzer, t, f¨
ur den Parameter, θ, sehr viele
∗
weiter Bootstrap-Sch¨atzer, t generiert. Die Anzahl der
Ziehungen bezeichnet man meist als R.
Aus diesen kann dann ein Vertrauensintervall berechnet
werden. Dazu gibt es verschiedene M¨oglichkeiten. Das sog.
“basic” Vertrauensintervall lautet
∗
∗
2t − t((R+1)(1−α/2))
, 2t − t((R+1)(α/2))
.
Parametrisches Bootstrapping
Aus der Stichprobe wird der Mittelwert x¯ und die Varianz σ 2
der Population gesch¨atzt.
Mit diesen Angaben wird die Verteilung beschrieben N (¯
x , σ2)
aus der dann wiederholt gezogen wird.
Das Parametrische Bootstrapping liefert oft bessere Resultate als
andere Verfahren.
Nichtparametrisches Bootstrapping
Das wiederholte Ziehen erfolgt direkt aus der Stichprobe.
Theoriefragen zu den bisherigen Inhalten
Was misst die Pearson bzw. die Spearman Korrelation?
Wann verwendet man eine Korrelation, wann eine Regression?
D¨
urfen die Messungen von verschiedenen Individuen in einer
Regression korreliert sein? Begr¨
unden Sie?
Warum transformieren wir Daten?
Theoriefragen zu den bisherigen Inhalten
Was bedeutet Overdispersal?
Wieso tritt bei der gew¨ohnlichen Regression nicht auch
Overdispersal auf?
Was untersucht der χ2 -Test oder der Fisher-Exakt-Test?
Welches sind die drei h¨aufigsten Arten wie man den
Zusammenhang zwischen Proportionen beschreibt?
Was ist ein Odds Ratio?
Theoriefragen zu den bisherigen Inhalten
Nennen Sie die drei Aspekte die man verwenden kann um ein
GLM zu charakterisieren.
Wie lautet die Linkfunktion einer logistischen Regression?
Wie gross ist der “dispersion parameter” in einem GLM?
Affiliation and copyright
Vielen Dank an Dr. Thomas Zumbrunn und Dr. Stefanie von Felten
von denen ich mehrere Ideen und Beispiele u
¨bernommen habe.
© Thomas Fabbro 2012
This work is licensed under the Creative Commons
Attribution-Noncommercial-Share Alike 2.5 Switzerland Licence.
To view a copy of this licence, visit
http://creativecommons.org/licenses/by-nc-sa/2.5/ch/
Document
Kategorie
Gesundheitswesen
Seitenansichten
16
Dateigröße
1 249 KB
Tags
1/--Seiten
melden