Statistik - uni-frankfurt.de · 1 Beschreibende Statistik Ein wesentlicher Teil der Statistik besch...

Statistik

Notizen zur Vorlesung im WS 2009/10

Gaby [email protected]

http://ismi.math.uni-frankfurt.de/schneider

Inhaltsverzeichnis

1 Beschreibende Statistik 11.1 Erwartungswert und zentraler Grenzwertsatz . . . . . . . . . . . . . . . 11.2 Quantile und Verteilungsfunktion . . . . . . . . . . . . . . . . . . . . . . 4

1.2.1 Ein allgemeingultiges Konfidenzintervall fur den Median . . . . . 61.2.2 Asymptotische Normalitat des Medians . . . . . . . . . . . . . . 71.2.3 Einschub: Lage- und Skalenparameter . . . . . . . . . . . . . . . 8

1.3 Ein Konfidenzband fur die Verteilungsfunktion . . . . . . . . . . . . . . 91.4 Graphische Darstellungen fur Verteilungen . . . . . . . . . . . . . . . . . 10

1.4.1 Das Histogramm . . . . . . . . . . . . . . . . . . . . . . . . . . . 101.4.2 Der Box-and-Whisker-plot (oder kurz: Boxplot) . . . . . . . . . . 101.4.3 Quantilplots . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13

2 Erste statistische Tests 142.1 Ein Permutationstest . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 142.2 Der Rangsummentest von Wilcoxon / Der Mann-Whitney-U-Test . . . . 142.3 Der z-Test . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 172.4 p-Wert, Signifikanz, Testmacht und Alternativhypothese . . . . . . . . . 18

3 Der t-Test 193.1 Zur Geometrie von Strichprobenmittelwert und -varianz . . . . . . . . . 19

3.1.1 Schatzer fur σ . . . . . . . . . . . . . . . . . . . . . . . . . . . . 193.1.2 Geometrie des Mittelwertes . . . . . . . . . . . . . . . . . . . . . 193.1.3 Geometrie der geschatzten Standardabweichung . . . . . . . . . . 21

3.2 Die t-Statistik . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 213.2.1 Die t-Verteilung und ihre Verwandten . . . . . . . . . . . . . . . 21

3.3 Die Standardnormalverteilung auf Rn und der Satz von Fisher . . . . . 233.4 Anwendung: Eine Stichprobe . . . . . . . . . . . . . . . . . . . . . . . . 243.5 Der zwei-Stichproben t-Test (fur gleiche Varianzen) . . . . . . . . . . . . 253.6 Zur Geometrie in der zwei-Stichproben-Situation . . . . . . . . . . . . . 27

4 Tests fur gepaarte Stichproben 284.1 t-Test fur gepaarte Stichproben . . . . . . . . . . . . . . . . . . . . . . . 28

4.1.1 Voraussetzungen und Spezialfall . . . . . . . . . . . . . . . . . . 294.2 Wilcoxon-Test fur gepaarte Stichproben/ Vorzeichenrangtest . . . . . . 304.3 Vorzeichentest . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 304.4 Permutationstest von Fisher . . . . . . . . . . . . . . . . . . . . . . . . . 31

5 Vergleich von mehr als 2 Stichproben und das normale lineare Modell 315.1 Vom t-Test zur einfaktoriellen Varianzanalyse . . . . . . . . . . . . . . . 31

5.1.1 Geometrie beim t-Test . . . . . . . . . . . . . . . . . . . . . . . . 315.1.2 Verallgemeinerung auf k(> 2) Gruppen . . . . . . . . . . . . . . 32

5.2 Das normale lineare Modell . . . . . . . . . . . . . . . . . . . . . . . . . 355.2.1 Testen linearer Hypothesen . . . . . . . . . . . . . . . . . . . . . 36

5.3 Paarweise Vergleiche, multiples Testen . . . . . . . . . . . . . . . . . . . 365.3.1 Bonferroni-Korrektur . . . . . . . . . . . . . . . . . . . . . . . . 375.3.2 Simultane Konfidenzintervalle nach Tukey . . . . . . . . . . . . . 38

5.4 Lineare Kontraste . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 395.5 Kruskal-Wallis-(H)-Test . . . . . . . . . . . . . . . . . . . . . . . . . . . 40

6 Lineare Regression 426.1 Modell und Parameterschatzer . . . . . . . . . . . . . . . . . . . . . . . 426.2 Die Geometrische Bedeutung von r und der Multiple Korrelationskoeffizient 44

6.2.1 Ein Konfidenzintervall fur die Korrelation . . . . . . . . . . . . . 456.3 Ein Konfidenzbereich fur β . . . . . . . . . . . . . . . . . . . . . . . . . 476.4 Multivariate Regression . . . . . . . . . . . . . . . . . . . . . . . . . . . 48

6.4.1 Orthogonale Regressoren . . . . . . . . . . . . . . . . . . . . . . 48

7 Schatzverfahren 497.1 Der Mittlere quadratische Fehler . . . . . . . . . . . . . . . . . . . . . . 507.2 Das Maximum-Likelihood Prinzip . . . . . . . . . . . . . . . . . . . . . . 52

7.2.1 Eigenschaften von Maximum-Likelihood-Schatzern . . . . . . . . 527.3 Die Momentenmethode . . . . . . . . . . . . . . . . . . . . . . . . . . . . 557.4 Kleinste Quadrate . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 567.5 Bootstrap-Konfidenzintervalle . . . . . . . . . . . . . . . . . . . . . . . . 567.6 Suffizienz und die Minimierung d. mittleren quadratischen Fehlers . . . 57

8 Likelihood-Quotienten-Tests 618.1 Einfache Hypothese und einfache Alternative . . . . . . . . . . . . . . . 628.2 Generalisierte Likelihood-Quotiententests . . . . . . . . . . . . . . . . . 63

9 Der χ2-Test 659.1 Normalverteilungen auf Rk . . . . . . . . . . . . . . . . . . . . . . . . . 679.2 Der χ2-Test . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68

1 Beschreibende Statistik

Ein wesentlicher Teil der Statistik beschaftigt sich mit Datenanalyse mit Hilfe stochas-tischer Modelle. Daten liegen oft vor in Form eines Zahlenvektors X = (x1, . . . , xn). Umdie Variabilitat dieser Beobachtungen zu modellieren, betrachtet man sie im stochasti-schen Modell als Realisierungen von ZVen X1, . . . , Xn. Wir schreiben daher wahlweiseX = (x1, . . . , xn) fur die Beobachtungen bzw. X = (X1, . . . , Xn) fur einen Vektor vonZufallsvariablen.

1.1 Erwartungswert und zentraler Grenzwertsatz

Definition Der Erwartungswert einer diskreten Zufallsvariablen X mit Werten in Sist definiert als

µ = EX =∑x∈S

xWs(X = x).

Der Erwartungswert einer Zufallsvariablen mit Wahrscheinlichkeitsdichte f(x) ist analog

µ = EX =∫R

xf(x)dx.

Schatzung des Erwartungswertes durch den Mittelwert Seien X1, . . . , Xn un-abhangige, identisch verteilte ZVe. Ein beliebter Schatzer fur den Erwartungswert istder Stichprobenmittelwert:

µ := X :=1n

n∑i=1

Xi.

Beobachtung 1: Erwartungstreue Fur X1, . . . , Xn u.i.v. mit E|X| < ∞ ist derMittelwert ein erwartungstreuer Schatzer fur den Erwartungswert:

E [µ] =1n

n∑i=1

EXi = EX1 = µ.

Beobachtung 2: Der mittlere quadratische Fehler Wie genau konnen wir µschatzen? Ein haufig verwendetes Kriterium fur die Gute eines Schatzers ist der mittlerequadratische Fehler: die erwartete quadratische Abweichung des Schatzers von dem zuschatzenden Wert. Schatzt man µ durch X, so ist der mittlere quadratische Fehler gerade

E((X − µ)2)EX=µ

= Var(X) = (1/n)Var(X).

Eine wichtige Unterscheidung: Die Standardabweichung der ZVe X ist√

Var(X), dieStandardabweichung des Mittelwerts von n unabhangigen ZVen Xi ist

Der Standardfehler : sd

(1n

∑i

Xi

)=

√Var(X)

n=

sd(X)√n

.

1

Beobachtung 3: Das Starke Gesetz der Großen Zahlen Sei X1, X2, . . . eine Folgeunabhangiger, identisch verteilter ZVen mit Erwartungswert µ. Xn := (1/n)

∑ni=1Xi.

Dann konvergiert Xn fast sicher gegen den Erwartungswert, d.h.

Ws(Xn → µ) = 1.

Das Starke Gesetz der Großen Zahlen impliziert stochastische Konvergenz des Stich-probenmittelwertes gegen den Erwartungswert: Fur große n liegt X mit hoher Wahr-scheinlichkeit nah beim Erwartungswert. Ein Schatzer, der stochastisch gegen den zuschatzenden Wert konvergiert, heißt konsistent. Der Stichprobenmittelwert ist also einkonsistenter Schatzer fur µ.

Die Verteilung von X Zwei Beispiele:

1. X1, . . . , Xn unabhangig und N (µ, σ2)-verteilt, also gemaß der Dichte

f(x) =1√

2πσ2e−

(x−µ)2

2σ2 , x ∈ R.

Dann ist der Stichprobenmittelwert normalverteilt mit Erwartungswert µ und Va-rianz σ2/n.

2. X1, . . . , Xn unabhangig und Exp(λ)-verteilt, also gemaß

f(x) = λe−λx, x > 0.

Der Stichprobenmittelwert hat Erwartungswert 1/λ und Varianz 1/(nλ2). Aber erist (zumindest fur kleine n) nicht normalverteilt. (In der Tat folgt die Summe von nunabhangigen Exp(λ)-verteilten ZVen einer Gammaverteilung mit den Parameternλ und n.)

Der Zentrale Grenzwertsatz Sei X1, X2, . . . eine Folge unabhangiger, identisch ver-teilter ZVen mit Erwartungswert µ und Varianz 0 < σ2 <∞. Sn :=

∑ni=1Xi. Dann ist

Sn approximativ normalverteilt, bzw.

Sn − nµσ√n−→ N (0, 1) fur n→∞ in Verteilung.

Ein asymptotisches Konfidenzintervall fur µ Der Zentrale Grenzwertsatz sagt,dass fur jede beliebige zugrundeliegende Verteilung der Werte X1, . . . , Xn die Verteilungdes Stichprobenmittelwertes (fur große n) einer Normalverteilung folgt. Damit konnenwir leicht ein erstes Konfidenzintervall konstruieren:Seien wiederum X1, . . . , Xn u.i.v. mit EX = µ und Var(X) = σ2 ∈ (0,∞). Sei µ =1/n

∑Xi. µ ist eine Zufallsvariable und schwankt um µ. In der Regel ist aber naturlich

µ 6= µ (und unsere Schatzung damit zwar gut, aber falsch). Wie falsch ist also unserSchatzer?Gesucht ist ein Intervall I(X), das mit hoher Wahrscheinlichkeit den wahren Erwar-tungswert uberdeckt, z.B.

Ws(I(X) 3 EX) ≥ 0.95. (1)

2

Ein solches Intervall nennt man ein 95%-Konfidenzintervall fur den Erwartungswert.Man beachte die Schreibweise: Nicht etwa der Erwartungswert ist hier zufallig, son-dern das Intervall! Nach dem ZGWS gilt (fur genugend großes n), dass µ approximativnormalverteilt ist. Das heißt, es gibt eine Zahl a so, dass

Ws(µ− aσ/√n ≤ µ ≤ µ+ aσ/

√n) = 0.95.

Bei der Normalverteilung mit Mittelwert µ und Varianz σ2 liegen etwa 95% der Masseim Intervall [µ− 2σ, µ+ 2σ] (Figur 1, der genaue Wert ist a = 1.96.) Also:

Ws(µ− 1.96σ/√

n ≤ µ ≤ µ+ 1.96σ/√

n) = 0.95.

Um dieses Intervall in ein Intervall fur µ zu transformieren, beobachten wir:

Ws(µ− 1.96σ/√n ≤ µ ≤ µ+ 1.96σ/

√n)

=Ws(−1.96σ/√

n ≤ µ− µ ≤ 1.96σ/√

n)=Ws(µ+ 1.96σ/

√n ≥ µ ≥ µ− 1.96σ/

√n)

Also erfullt I(X) = [µ− 1.96σ/√n, µ+ 1.96σ/

√n] Bedingung (1).

Normalverteilung

µµ

Dic

hte

µµ ++ σσµµ −− σσ

68%

µµ ++ 2σσµµ −− 2σσ

95%

µµ ++ 3σσµµ −− 3σσ

99.7%

Abbildung 1: Normalverteilung: Flache in Intervallen [µ− kσ, µ+ kσ], k = 1, 2, 3.

Was kann schief gehen?

• Die Voraussetzungen konnten nicht zutreffen, z.B. existiert EX nicht fur alle Ver-teilungen. Wir betrachten die Standard-Cauchy-Verteilung, mit Dichte

f(x) =1π

11 + x2

, x ∈ R.

Hier ist zwar f(x) integrierbar, aber nicht xf(x), so dass der Erwartungswert nichtexistiert. Das starke Gesetz der Großen Zahlen gilt hier nicht — die normiertenPartialsummen Xn = 1/n

∑ni=1Xi konvergieren nicht (Figur 2).

3

−4 −2 0 2 4

0.0

0.2

0.4

x

Dic

hte

Normalvert.Cauchyvert.

0 2000 4000 6000 8000 10000

−4

−2

02

n

Sn

n

Abbildung 2: Fehlende Konvergenz des Mittelwertes bei der Cauchy-Verteilung.

• Die Aussage des ZGWS ist asymptotisch, aber die Gute der Asymptotik hangt vonder zugrundeliegenden Verteilung ab. Fur manche Verteilungen muss fur approxi-mative Normalitat n erst sehr groß werden. Zum Beispiel, wenn eine Verteilungwinziges Gewicht auf sehr großen Werten hat:X habe Verteilung ν0 mit Erwartungswert µ.

Xε :=

{X mit W’keit 1− ε1/ε2 mit W’keit ε

µε := EXε = (1− ε)µ+ ε1ε2−→∞ fur ε→ 0.

Fur festes n gilt fur jedes von X = (X1, . . . , Xn) abhangige Intervall I(X):

Ws(I(X1, . . . , Xn) = I(Xε1 , . . . , X

εn)) ε→0−→ 1.

Konfidenzintervalle sollten moglichst klein sein, insbesondere beschrankt. Daherfolgt mit

Ws(I(X1, . . . , Xn) 3 µε)ε→0−→ 0

Ws(I(Xε1 , . . . , X

εn) 3 µε)

ε→0−→ 0.

Die Approximation aus dem ZGWS kann also beliebig schlecht werden.

Fur so pathologische Falle, z.B. bei ausreißerbehafteten Daten, oder wenn man wenigeInformationen uber die zugrundeliegende Verteilung hat, verwendet man gerne nichtpa-rametrische Ansatze, z.B. Quantile.

1.2 Quantile und Verteilungsfunktion

Definition Die V erteilungsfunktion einer ZVen X ist definiert als

F (x) = Ws(X ≤ x) = ν((−∞, x]),

wobei ν die Verteilung von X bezeichnet.Diskret: F (x) =

∑y≤x Ws(X = y).

Hat X Dichte f(x): F (x) =x∫−∞

f(y)dy.

4

Schatzung der Verteilungsfunktion Zu einem Datenvektor X = (x1, . . . , xn) istdie empirische V erteilungsfunktion definiert als

Fn(X) =1n|{xi|xi ≤ x}|.

X=(1,2,4,7,9)

x

#xi≤≤

x

0.0

0.4

0.8

−2 0 2 4 6 8 10 12

●

●

●

●

●

●

●

●

●

●

100 normalverteilte ZVen

x

#xi≤≤

x

0.0

0.4

0.8

0 2 4 6 8 10

Abbildung 3: Beispiele fur empirische Verteilungsfunktionen

Definition Fur 0 < p < 1 heißt eine reelle Zahl qp p−Quantil der Verteilung ν, wenn

Ws(X ≤ qp) ≥ p und Ws(X ≥ qp) ≥ 1− p.

Bsp N (0, 1) (Fig. 4 A) Hat F Dichtefunktion f = F ′, die im Bereich {x ∈ R|0 <F (x) < 1} strikt positiv ist (d.h. F ist streng monoton steigend), so ist ∀p ∈ (0, 1) dasp-Quantil eindeutig bestimmt.

Bsp Bi(2,1/2) (Fig. 4 B) q1/2 = 1, q1/6 = 0, q1/4 =?Alle Werte im Intervall [0, 1] sind 25%-Quantile der Verteilung Bi(2,1/2). Also:

Die Menge aller p-Quantile einer Verteilung ν ist ein Intervall [qp,1, qp,2], wobei

qp,1 := sup{x ∈ R|F (x) < p} und qp,2 := inf{x ∈ R|F (x) > p}.

N(0,1)

x

P(X

≤x)

0.0

0.4

0.8

−3 −2 −1 0 1 2 3

Bi(2,1/2)

x

P(X

≤x)

0.0

0.4

0.8

−1 0 1 2 3

●

●

●

●

●

●

Abbildung 4: Verteilungsfunktionen von N (0, 1) und Bi(2,1/2)

5

Spezielle Quantile

q0.5 ein Median (’der Median’:q0.5,1 + q0.5,2

2)

q0.25 ein 1. Quartil (’das 1. Quartil’:q0.25,1 + q0.25,2

2)

q0.75 ein 3. Quartil (’das 3. Quartil’:q0.75,1 + q0.75,2

2)

Schatzer fur den Median Sei (X(1), X(2), . . . , X(n)) mit X(1) ≤ X(2) ≤ . . . ≤ X(n)

die Ordnungsstatistik, d.h. die der Große nach geordnete Stichprobe. Dann wird derMedian geschatzt durch den Stichprobenmedian

m =

{X(k+1) fur n = 2k + 11/2(X(k) +X(k+1)) fur n = 2k.

1.2.1 Ein allgemeingultiges Konfidenzintervall fur den Median

Im Gegensatz zum Mittelwert kann man fur den Median unabhangig von der zugrunde-liegenden Verteilung sogar fur festes n ein Intervall I(X) angeben so, dass

Ws(I(X) 3 m) moglichst groß (z.B. 0.95).

Wir verwenden dazu wieder die Ordnungsstatistik (X(1), . . . , X(n)) und beginnen miteinem Beispiel:

Beispiel: Eine Stichprobe der Große n = 11. Betrachte zunachst das Intervall I1 :=[X(1), X(n)]. Es gilt:

Ws(m /∈ I1) = Ws(X(1) > m) + Ws(X(n) < m)

= 2 ·(

12

)11

≈ 0.001.

Konnen wir dieses Intervall verkleinern so, dass Ws(m /∈ I) gerade noch kleiner ist als0.05? Was ist mit I2 = [X(2), X(n−1)]?

Ws(m /∈ I2) = Ws(X(2) > m) + Ws(X(n−1) < m)

= Ws(hochstens ein Xi ist kleiner m) + Ws(hochstens ein Xi ist großer m)

= 2

[(111

)(12

)11

+(

110

)(12

)11]≈ 0.012

in R: = 2 ∗ pbinom(1, 11, 0.5)

Um ein moglichst kleines Konfidenzintervall zu finden so, dass Ws(m ∈ I) ≥ 0.95, suchenwir eine Zahl ` so, dass 2pbinom(`, 11, 0.5) gerade noch kleiner ist als 0.05. Das gehtz.B. in R: mit der implementierten Quantilfunktion: ` = qbinom(0.05/2, 11, 0.5). Die

6

Wahrscheinlichkeit 0.05 nennt man das Niveau des Konfidenzintervalls. Allgemein alsofur n und Niveau α

` = qbinom(α/2, n, 0.5).

Ahnliche Uberlegungen erlauben die Konstruktion von Konfidenzintervallen fur jedesbeliebige Quantil und jede zugrundeliegende Verteilung fur festes n. (Ubung)

1.2.2 Asymptotische Normalitat des Medians

Fur symmetrische Verteilungen fallen Erwartungswert und Median zusammen. Was soll-te man also nehmen, um die Lage zu schatzen? Nach dem ZGWS ist der Stichproben-mittelwert asymptotisch normalverteilt mit Varianz σ2/n. Fur den Median gibt es eineahnliche Aussage:

Satz Seien X1, . . . , Xn u.i.v. mit einer differenzierbaren Verteilungsfunktion F . Sei mder Median der Verteilung, und es gelte F ′(m) > 0 (in diesem Fall ist der Medianeindeutig). Sei Mn der Stichprobenmedian. Dann gilt

√n(Mn −m) −→ N (0, (2F ′(m))−2) in Verteilung.

Das heißt, der Stichprobenmedian ist approximativ normalverteilt, und die approxima-tive Varianz ist umso kleiner, je großer der Wert der Dichte am Median ist.Beweis: Es ist zu zeigen, dass

Ws(√n(Mn −m) ≤ a) −→Ws(Z ≥ −2aF ′(m)) fur n→∞,

fur eine Standardnormalverteilte Zufallsvariable Z. Wir definieren dazu Indikatorvaria-blen

Yi :=

{1 falls Xi ≤ m+ a/

√n

0 sonst

Damit giltWs(Yi = 1) = F (m+ a/

√n) =: pn.

Wir konnen Ws(Mn ≤ m+ a√n) nun zuruckfuhren auf die binomialverteilte ZVe

∑Yi,

die approximativ normalverteilt ist:

Ws(Mn ≤ m+ a/√n)

n ungerade= Ws

(∑Yi ≥

n+ 12

)= Ws

( ∑Yi − npn√

npn(1− pn)≥ (n+ 1)/2− npn√

npn(1− pn)

)

7

Nach dem ZGWS konvergiert Wn :=PYi−npn√

npn(1−pn)in Verteilung gegen eine standardnor-

malverteilte ZVe Z. Bei Grenzwertbetrachtung des rechten Terms, cn, erhalten wir

cn :=(n+ 1)/2− npn√

npn(1− pn)=−n(pn − 1/2)√npn(1− pn)

+1/2√

npn(1− pn)

=1√

pn(1− pn)−(pn − 1/2)

1/√n

+1/2√

npn(1− pn)

=1√

pn(1− pn)−(F (m+ a/

√n)− F (m))

1/√n

+1/2√

npn(1− pn)

=−a√

pn(1− pn)︸︷︷︸→−2a

F (m+ a/√n)− F (m)

a/√n︸︷︷︸

→F ′(m)

+1/2√

npn(1− pn)︸︷︷︸→0

n→∞−→ −2aF ′(m)

Damit sehen wir mit c = −2aF ′(m)

|Ws(Wn ≤ cn)−Ws(Z ≤ c)| ≤ |Ws(Wn ≤ cn)−Ws(Z ≤ cn)|+|Ws(Z ≤ cn)−Ws(Z ≤ c)|

Der zweite Summand konvergiert gegen 0, weil cn → c gilt und die VerteilungsfunktionFZ von Z stetig ist. Um zu sehen, dass der erste Summand gegen 0 konvergiert, erinnernwir uns, dass aus Konvergenz in Verteilung fur stetige Verteilungsfunktion FZ auchgleichmaßige Konvergenz folgt. Weil also FWn gleichmaßig gegen FZ konvergiert, gilt|Ws(Wn ≤ cn)−Ws(Z ≤ cn)| ≤ supx∈R |Ws(Wn ≤ x)−Ws(Z ≤ x)| → 0. Also gilt

Ws(Wn ≤ cn)→Ws(Z ≤ c).

Damit erhalten wir wie gewunscht

Ws(√n(Mn −m) ≤ a) −→Ws(Z ≥ −2aF ′(m))

Wegen Ws(Z ≥ −2aF ′(m)) = Ws(Z/2F ′(m) ≤ a) konvergiert daher Mn gegen einenormalverteilte ZVe mit Varianz 1/(2F ′(m))2. �

1.2.3 Einschub: Lage- und Skalenparameter

Lageparameter Mittelwert, Median und die Quantile im allgemeinen heißen Lagepa-rameter, weil sie bei einer Verschiebung der Verteilung um θ ebenfalls um θ verschobenwerden:Seien µ der Mittelwert und qp das p-Quantil der Verteilung ν mit Dichte f0(x) undVerteilungsfunktion F0(x), so ist µ + θ der Mittelwert und qp + θ das p-Quantil derVerteilung νθ mit Dichte fθ(x) = f0(x− θ) und Verteilungsfunktion Fθ(x) := F0(x− θ).Parameter mit dieser Eigenschaft heißen Lageparameter einer Verteilung.

Skalenparameter beschreiben die Breite einer Verteilung. Ein Parameter θ heißt Ska-lenparameter einer Verteilung ν, wenn gilt:Hat X die Verteilung ν(0,1) mit Dichte f(0,1)(x) und Verteilungsfunktion F(0,1)(x), so

8

hat Y := θX die Verteilung ν(0,θ) mit Verteilungsfunktion F(0,θ)(x) := F(0,1)(x/θ) undDichte f(0,θ)(x) = 1/θf(0,1)(x/θ).Beispiel: Die Standardabweichung σ =

√E(X − µ)2 ist ein Skalenparameter, denn hat

X Standardabweichung σ, so hat Y := cX Standardabweichung cσ.

1.3 Ein Konfidenzband fur die Verteilungsfunktion

Definition Die Quantiltransformation. Fur 0 < p < 1 sei

F−1(p) := min{x ∈ R|F (x) ≥ p} = sup{x ∈ R|F (x) < p} = qp,1

dieQuantilfunktion (verallgemeinerte Umkehrfunktion von F , das minimale p-Quantil).

Bsp. 1 Bei Verteilungen mit endlichem Trager reduziert sich F−1(p) jeweils auf dieuntere Intervallgrenze aller p-Quantile.

Bsp. 2 Bei Verteilungen mit positiver Dichte f = F ′, z.B. Exponentialverteilung:

Fλ(x) =

{1− e−λx x ≥ 00 sonst

F ist fur x ≥ 0, also fur 0 < p < 1 umkehrbar, also gilt

F−1(p) = − 1λ

log(1− p) fur p ∈ (0, 1).

Lemmaa) U sei uniform verteilt auf [0, 1], d.h. Ws(U ∈ [a, b]) = b− a fur 0 ≤ a ≤ b ≤ 1. Danndefiniert

X := F−1(U)

eine ZVe mit Verteilungsfunktion F .Bew.:Nach der Definition von F−1(u) := min{x ∈ R|F (x) ≥ u} gilt wg Ws(U ∈ (0, 1)) = 1

F−1(U) ≤ r ⇐⇒ F (r) ≥ U, alsoWs(X ≤ r) = Ws(U ≤ F (r)) = F (r).

b) U1, . . . , Un seien unabhangig und uniform auf [0, 1]. Es seien

Gn(v) :=1n|{i|Ui ≤ v}| die empir. Verteilungsfunktion der Ui,

Fn(r) :=1n|{i|Xi ≤ r}| die empir. Verteilungsfunktion der Xi.

Dann ist ∀r Gn(F (r)) = 1n |{i|Ui ≤ F (r)}| = 1

n |{i|F−1(Ui) ≤ r}| = Fn(r), also ist

(Gn(F (r))) genauso verteilt wie (Fn(r)).

Korollar X sei ZVe mit stetiger Verteilungsfunktion F , dann folgt F (X) ∼ Unif[0, 1].

9

Korollar (Kolmogorov-Smirnov)Sei Gn(v) := 1

n |{i|Ui ≤ v}| wie gehabt und c ≥ 0. Dann gilt:

Ws(supr∈R|Fn(r)− F (r)| ≤ c) = Ws( sup

u∈F (r)|Gn(u)− u| ≤ c) ≥Ws( sup

u∈[0,1]|Gn(u)− u| ≤ c)

Sei alsocn(α) := min{c ≥ 0|Ws( sup

u∈[0,1]|Gn(u)− u| ≤ c) ≥ 1− α}.

Dann ist mitWs(F (r) ∈ [Fn(r)± cn(α)]∀r) ≥ 1− α

ein (1− α)-Konfidenzband fur den Graphen von F gegeben.

Werte fur cn(α) kann man Tabellen entnehmen. Es ist bekannt, dass

cn(α) ≤√− log(α/2)

2n.

1.4 Graphische Darstellungen fur Verteilungen

Wir wollen kurz ein paar verbreitete Darstellungsformen fur Daten und ihre Implemen-tierung in R kennen lernen (vgl. Fig. 5).

1.4.1 Das Histogramm

Ein Histogramm liefert eine umfangreiche Darstellung einer Haufigkeitsverteilung: Wieviele Beobachtungen fallen in welches Intervall? Ein Histogramm, auf Flache 1 normiert,stellt eine rudimentare Schatzung fur die Dichte dar. Die Hohe eines Balkens gibt denAnteil der Verteilung pro Einheit der x-Achse wieder (und ist deswegen von der Wahlder Einheit abhangig).

1.4.2 Der Box-and-Whisker-plot (oder kurz: Boxplot)

Einen Zahlenvektor x1, . . . , xn kann man kurz in funf Zahlen beschreiben: Minimum,25%-Quantil (Q1, das erste Quartil), 50%-Quantil (Median), 75%-Quantil (Q3, dasdritte Quartil), und Maximum. Ein Boxplot liefert eine graphische Darstellung einerahnliche Kurzubersicht: 25%-und 75%-Quantil bilden die Rander einer Box um den Me-dian. Hier liegen also die mittleren 50% der Beobachtungen. Der InterquartilabstandQ3-Q1 ist neben der Standardabweichung eine gangige Art, die Breite einer Verteilungzu beschreiben.Die Whisker (‘Barthaare’) konnen bis zu den Extremwerten reichen. Aber in der Regelbetragt ihre Lange das 1.5-Fache des Interquartilabstands Q3-Q1, wobei sie nur bis zurextremsten Beobachtung innerhalb dieses Bereiches reichen. Werte außerhalb werdendurch Punkte dargestellt.

10

Tabelle 1: Histogramme in R

Gegeben sei ein Datenvektor x<-c(4,2,5.5,0.7,...)Histogrammbefehl: hist(x)nutzliche OptionenIntervalle: hist(x,breaks=12) produziert 12 Balken

hist(x,breaks=seq(min(x),max(x),by=0.1)) Balken der Breite 0.1

Normiert die Flache auf 1: hist(x,freq=FALSE)speichere Ausgabe h<-hist(x,plot=FALSE)weitere Informationen mit ?hist

Auch fur andere Darstellun- gen gultige Optionen:Titel hist(x,main=‘‘Histogramm des Datenvektors x’’)Achsenbeschriftung hist(x,xlab=‘‘Große’’, ylab=’’Anzahl’’)dargestellter Bereich hist(x,xlim=c(0,3), ylim=...)Farbe: hist(x,col=’’yellow’’) (Alle Farbnamen: colors())

Tabelle 2: Boxplots und stripcharts in R

Fur einen Datenvektor x<-c(4,2,5.5,0.7,...)einfach boxplot(x)

fur zwei Boxplots nebeneinander:Ein Datenvektor x<-c(4,2,5.5,0.7,...)Ein Vektor, der die Gruppen definiert Gruppe<-c(0,0,0,...,1,1,1,...)und dann einfach boxplot(x∼ Gruppe)

ganz analog fur stripchart: stripchart(x∼ Gruppe)nutzliche Optionen fur stripchart:‘ruttle’ die Datenpunkte ein wenig: stripchart(x∼ Gruppe,method=’jitter’)Beschriftung der Gruppen ...,group.names=c(’Gr. 1’,’Gr. 2’))

11

Tabelle 3: Graphikparameter in R

par(...) setzt Graphikparametercex.axis Große von Achsentitelncex.lab Große von Achsenbeschriftungencex.main Große der Uberschriftmfcol=c(2,3) Grafikfenster mit 2 Zeilen, 3 Spalten?par weitere nutzliche Parameter

1000 std.normalverteilte ZVen

x

Anz

ahl

−3 −2 −1 0 1 2 3

020

4060

80

●

Gruppe 1 Gruppe 2

−2

02

46

−2 0 2 4 6G

rupp

e 1

Gru

ppe

2

●

●● ●

●

●●

● ●

●

●

●

●

●

●

●

●

●

●●

●

●

●

●

●

●

●

●

●

●

●

●

●

●●

●

●

●

●●

●●

●

●●

●

● ●

●

●

●●●●

●

●

●

●● ●

●●

●●● ●

●

●

●

●

●

●

●

●●

●

●

●●

●●

●

●

●

●

●●

●

●

●

●●

●

●

●●

●● ●●

●

●

●

●

●

● ●●

●

●

●

●

●● ●

●

●

●●● ●

●

●

●

●

●

●

●

●

●

●

●●

●

●

● ●

●●

●

● ●

●

●

●

●

●

●

●

●

●

●

●

●●●

●

●●

●●●

●

●

●●

●

●

●●

●

●

●

●●

●

●

●

● ●

●

●●

●

●

●

●

●

●

●

●

●

● ●

●

●

●

●●

●

Abbildung 5: Beispiele fur Histogramm, Boxplot und Stripchart# Graphikparameterpar(mfrow=c(1,3),cex.lab=1.7,cex.axis=1.4,cex.main=1.8,font.main=1,mar=c(4,4.8,3,1))# Graphik 1hist(rnorm(1000),col=’yellow’,breaks=30,xlab=’x’,ylab=’Anzahl’,main=’1000std.normalverteilte ZVen’)# Datenerzeugung fur Graphiken 2 und 3x<-rep(c(’Gruppe 1’,’Gruppe 2’),each=100)y<-c(rnorm(100),rnorm(100,mean=2,sd=2))# Graphik 2boxplot(y∼x,col=’red’)# Graphik 3stripchart(y x,col=’blue’,pch=19,method=’jitter’)

12

1.4.3 Quantilplots

Wie ist die Form einer empirischen Verteilung? Vielleicht ungefahr wie eine Normal-verteilung (mit entsprechendem Mittelwert und Varianz)? Oder hat sie z.B. schwerereSchwanze, ist sie etwa schief, oder hat sie Locher?

qqnorm Um eine empirische Verteilung mit einer Normalverteilung zu vergleichen,konnte man entweder eine Normalverteilung anpassen (d.h. Mittelwert und Varianzschatzen) und als Kurve uber ein Histogramm der empirischen Haufigkeitsverteilung le-gen. Oder man verwendet einen Normalplot : Man identifiziert die Zahlen x1, . . . , xn mitihren Quantilen (in R werden die Zwischenschritte ausgewertet: 1/(2n), 3/(2n), . . . , (2n−1)/(2n)) und tragt die entsprechenden Quantile der Normalverteilung gegen die Beob-achtungen auf (Figur 6, unten links).

qqplot Dasselbe kann man machen, um zwei empirische Verteilungen zu vergleichen:Man identifiziert die Beobachtungen mit ihren Quantilen und tragt sie gegeneinanderauf. Fur zwei Datensatze x1, . . . , xn und y1, . . . , yn mit demselben Stichprobenumfangtragt man jeweils die i-t-kleinsten Beobachtungen x(i), y(i) gegeneinander auf (Figur 6,unten rechts).

Dichte

−2 0 2

Verteilungsfunktion

0.0

0.5

1.0

−2 0 2

qqnorm

Beobachtungen: N(5,4)

Qua

ntile

d. N

orm

alve

rtei

lung

0 5 10

−2

02

●

●

●

● ● ●

●

●●● ●

●●

●●

●

●

● ●

●

qqplot

St.probe aus Exp(1)

St.p

robe

aus

N(0

,1)

0.0 0.5 1.0 1.5 2.0 2.5 3.0

−1.

00.

01.

02.

0

Abbildung 6: Quantilplots: Erstellung und Interpretation

13

2 Erste statistische Tests

2.1 Ein Permutationstest

Hat das Horen von (mutterlich klingenden) Herztonen einen Effekt auf die Gewichtszu-nahme von Neugeborenen?Wir betrachten eine Untersuchung von Lee Salk (The Role of the Heartbeat in the Re-lation between Mother and Infant, Scientific American, May 1973, p. 26-29). Aus einerZahl Neugeborener wurden zufallig einige ausgewahlt, die in ihren ersten Lebenstagenohne Unterbrechung den Ton des Herzschlags eines Erwachsenen horen durften. Sonstwurden alle Babies gleich behandelt. Beobachtet wurde die Gewichtszunahme vom ers-ten auf den vierten Lebenstag.Hat die Behandlung (Horen der Herztone) einen Effekt auf die Gewichtszunahme? Gibtes einen signifikanten Unterschied in der Gewichtszunahme zwischen behandelten undunbehandelten Babies? Oder konnten die beobachteten Unterschiede durch Zufall zu-stande gekommen sein, auch wenn die Behandlung keinen Effekt hat?Wir haben aus m+ n Individuen zufallig m gezogen und behandelt:

x1, . . . , xm in der Behandlungsgruppe, mit x =1m

∑xi,

y1, . . . , yn in der Kontrollgruppe, mit y =1n

∑yi.

Wir beobachten x > y. Wie groß muss der Unterschied sein, damit wir auch meinen,das gilt ebenso fur die Populationsparameter µx > µy? Wir betrachten die Nullhypo-these µx = µy, dass die Behandlung keinen Effekt hat. Unter der Nullhypothese sinddie x1, . . . , xm einfach eine zufallige Stichprobe aus den Geburtsgewichtszunahmen derm+n Babies, denn dann hatte die Aufteilung der m+n Zahlen in die beiden Gruppengenauso gut anders ausfallen konnen. Also konnen wir x vergleichen mit dem Wert einerzufallig aus den m + n Zahlen gezogenen (also einer durch Permutation gewonnenen)Stichprobe:Ziehe rein zufallig (ohne Zurucklegen) Z1, . . . , Zm aus der Liste x1, . . . , xm, y1, . . . , ynund berechne die Wahrscheinlichkeit Ws

(1m

∑Zi ≥ x

).

Einfacher: Approximiere diese Wahrscheinlichkeit durch Simulation: Ziehe (z.B.) 10000mal Z1, . . . , Zm und notiere die relative Haufigkeit fur das Ereignis

{1m

∑Zi ≥ x

}.

Falls in weniger als 5% der Falle 1m

∑Zi ≥ x, so kann man sagen: Die Unterschiede

zwischen den Gruppen waren auf dem 5%-Niveau signifikant. Bedeutung: Wenn es kei-nen Behandlungseffekt gabe, beobachten wir solche Unterschiede in weniger als 5% derVersuche durch Zufall.

2.2 Der Rangsummentest von Wilcoxon / Der Mann-Whitney-U-Test

Der Permutationstest ist fur viele Fragestellungen einsetzbar, erfordert nur die Wahleiner adaquaten Teststatistik (hier: 1/m

∑Zi) und die Simulation mit dem Rechner.

Der Test macht keine Annahmen uber die zugrundeliegenden Verteilungen und ist da-her allgemein verwendbar. Da aber die Verteilung der Testsatistik 1/m

∑Zi auch von

den beobachteten Werten x1, . . . , xm, y1, . . . , yn abhangt, ist sie nur durch Simulation zuermitteln.

14

Will man diese Abhangigkeit von den Werten loswerden, kann man z.B. die Rangsum-menstatistik von Wilcoxon verwenden. Idee: Statt der Werte werden ihre Range verwen-det:Beispiel 1: Wenn alle Werte unterschiedlich sind

Beobachtung x 2 0 2.9 2.5 5Rang von x 2 1 4 3 5

Beispiel 2: Bei Rangbindungen (engl.: ties)

Beobachtung x 2 0 2.9 2 5Rang von x 2.5 1 4 2.5 5

Vergleicht man (wie im Beispiel der geherzten Babies) die Beobachtungen aus zweiGruppen, so werden die Range uber alle Beobachtungen gebildet. Beispiel:

Gruppe 1 1 1 2 1 2 1 2 2 2Beobachtung x −2 0 0.5 0.7 1.3 1.8 2.8 4 4.2 4.3Rang von x 1 2 3 4 5 6 7 8 9 10

Rangsumme in Gruppe 1: 18, Rangsumme in Gruppe 2: 10(10+1)2 − 18.

Um die Rangsummen bewerten zu konnen, mussen wir wissen, wie sie unter der Null-hypothese des rein zufalligen Ziehens von m Zahlen Z1, . . . , Zm aus den n + m Zahlenx1, . . . , xm, y1, . . . , yn verteilt waren. In dem Fall ist

(Rang(Z1), . . . ,Rang(Zm)) d= (R1, . . . , Rm),

wobei Ri rein zufallige Zuge aus {1, . . . ,m+n} ohne Zurucklegen sind. Die Rangsummebezeichnen wir mit

S :=m∑i=1

Ri.

Wir betrachten im Folgenden den Fall ohne Rangbindungen, in dem alle n + m Werteverschieden sind.

Bemerkung 1 Unter der Nullhypothese rein zufalligen Ziehens ist

ES =m(m+ n+ 1)

2,

denn

ERi =m+n∑k=1

k · 1m+ n

=(m+ n)(m+ n+ 1)

2· 1m+ n

=m+ n+ 1

2,

ES = m · ER1 =m(m+ n+ 1)

2.

15

Bemerkung 2 Unter der Nullhypothese ist

Var(S) =mn(m+ n+ 1)

12.

Beweis:Wir betrachten wieder den Fall, in dem alle n+m Werte verschieden sind, und konnenuns daher auf eine angenehmere Darstellung von S zuruckziehen:Seien U1, . . . , Um, V1, . . . , Vn unabhangig und gleichverteilt auf [0, 1]. Die Rangsummeder Ui, S =

∑Rang(Ui), ist dann verteilt wie die Rangsumme S =

∑Ri von oben.

Es reicht also, die Behauptung fur S zu zeigen. Der Rang eines Wertes Ui setzt sichzusammen aus der Anzahl der Uj , die kleiner oder gleich Ui sind, und der Anzahl derVj , die kleiner Ui sind:

Rang(Ui) =m∑k=1

I{Uk≤Ui} +n∑j=1

I{Vj<Ui}

=⇒ S =m∑i=1

Rang(Ui) =m(m+ 1)

2+

m∑i=1

n∑j=1

I{Vj<Ui}.

Damit ist

Var(S) = Var

m∑i=1

n∑j=1

I{Vj<Ui}

=∑(i,j)

Var(I{Vj<Ui}

)+

∑(i1, j1), (i2, j2)i1 6= i2, j1 6= j2

Cov(I{Vj1<Ui1}, I{Vj2<Ui2}

)

+m∑i=1

∑(j1,j2)

Cov(I{Vj1<Ui}, I{Vj2<Ui}

)+

n∑j=1

∑(i1,i2)

Cov(I{Vj<Ui1}, I{Vj<Ui2}

)Der zweite Summand ist Null, weil fur i1 6= i2 und j1 6= j2 die Ereignisse {Vj1 < Ui1}und {Vj2 < Ui2} unabhangig sind. Fur den dritten und vierten Summanden gilt

Cov(I{Vj1<Ui}, I{Vj2<Ui}

)= E

[I{Vj1<Ui} · I{Vj2<Ui}

]− 1

4

= Ws({Vj1 < Ui} ∩ {Vj2 < Ui})−14

=13− 1

4=

112

Damit ist

Var(S) =nm

4+mn(n− 1)

12+mn(m− 1)

12=mn(m+ n+ 1)

12.

�

Die Rangsummenstatistik S ist eine Summe von ZVen, von denen fur große m,n dieuberwiegende Anzahl unabhangig sind (fur die meisten Paare (i1, j1), (i2, j2) gilt i1 6= i2und j1 6= j2 und damit sind {Vj1 < Ui1} und {Vj2 < Ui2} unabhangig). Dies macht dasfolgende Resultat einleuchtend (das aber hier nicht bewiesen wird):

16

Satz (Hoeffding):Fur m,n→∞ gilt:

Sn,m − ESn,m√Var(Sn,m)

→ N (0, 1) in Verteilung.

Der Rangsummentest von Wilcoxon in R Zum Vergleich von zwei Vektoren xund y verwendet man den Befehl wilcox.test(x,y). R berechnet die Mann-Whitney-Statistik

U =m∑i=1

n∑j=1

I{Vj<Ui},

die wegen S = m(m+1)2 + U aquivalent ist zur Rangsummenstatistik S. Ohne weitere

Spezifikation wird fur m + n ≤ 50 die exakte Verteilung von U ermittelt und damitder exakte p-Wert. Fur m + n > 50 wird die Normalapproximation verwendet (es seidenn, man spezifiziert wilcox.test(x,y,exact=TRUE). Fur diesen Fall findet man inder Hilfe (?wilcox.test) eine Warnung: This function can use large amounts of memoryand stack (and even crash R if the stack limit is exceeded) if exact = TRUE and onesample is large.)

Achtung Der Wilcoxon-Test testet die Nullhypothese, dass beide Verteilungen gleichsind, gegen die Alternative, dass die eine gegenuber der anderen verschoben ist. Erist nicht sensitiv fur andere Arten von Alternativen (z.B. wenn sich die Streuungenunterscheiden, nicht aber die Mittelwerte).

Die verbreitetste Fragestellung vergleicht die Erwartungswerte zweier Verteilungen.Sie ist am leichtesten mit Hilfe eines z-Tests zu behandeln, sofern die Standardabwei-chung bekannt ist:

2.3 Der z-Test

Aus dem asymptotischen Konfidenzintervall fur den Mittelwert, das wir in Abschnitt1.1 kennengelernt haben, konnen wir leicht einen ersten statistischen Test konstruieren.

Ein Beispiel Wir mochten untersuchen, ob sich das Geburtsgewicht von Babies ineiner bestimmten Region Deutschlands von dem Geburtsgewicht im gesamten Bundes-gebiet unterscheidet. Angenommen, wir wussten (z.B. vom statistischen Bundesamt),dass das mittlere Geburtsgewicht in Deutschland µ = 3.3 kg betragt und eine Stan-dardabweichung von σ = 0.5 kg hat. Um zu untersuchen, ob sich das Geburtsgewichtin der ausgewahlten Region davon unterscheidet, beobachten wir eine Stichprobe vonn = 50 Babies, die dort geboren wurden. Ihr mittleres Geburtsgewicht betragt x = 3.5kg. Wie stark weicht x von µ ab? Spricht diese Abweichung dafur, dass das mittlereGeburtsgewicht µR in dieser Region tatsachlich hoher ist, oder haben wir es mit einerzufalligen Schwankung zu tun?Nach dem ZGWS bewerten wir Abweichungen des Mittelwertes in Einheiten des Stan-dardfehlers, σ/

√n, unsere Teststatistik ist also

z =∣∣∣∣ x− µσ/√n

∣∣∣∣ .17

Falls das Geburtsgewicht in der betreffenden Region sich von der ubrigen Populationnicht unterscheidet (also unter der Nullhypothese µR = µ), ist z nach dem ZGWSapproximativ standardnormalverteilt. Auf dieser Basis beurteilen wir, wie wahrscheinlichdas beobachtete Ereignis unter der Hypothese µR = µ ist.

In unserem Beispiel ist z = 2.83. Fur eine N (0, 1)-verteilte ZVe Z ist aber

Ws(|Z| > 2.83) < 0.005.

Unter der Hypothese µR = µ haben wir also ein sehr unwahrscheinliches Ereignis beob-achtet. Man sagt auch: Die Daten sprechen (auf dem 0.5%-Niveau) gegen die HypotheseµR = µ.

2.4 p-Wert, Signifikanz, Testmacht und Alternativhypothese

Der p-Wert bezeichnet die Wahrscheinlichkeit, ein mindestens so extremes Ereigniszu beobachten, wenn (d.h. obwohl) die Nullhypothese zutrifft. Im Beispiel von eben alsop = 0.005. Man spricht ublicherweise von statistischer Signifikanz, wenn p < 0.05 (gerneauch durch * markiert), hochsignifikant bedeutet oft p < 0.01 (**) oder p < 0.001 (***).Nicht signifikant (‘n.s.’) heißt entsprechend: p ≥ 0.05.

Wie ist ein nicht signifikantes Ergebnis zu interpretieren? Wie sehr spricht es fur dieNullhypothese und gegen die Alternative? Das hangt von der gewahlten Alternativhy-pothese ab.

Beispiel X1, . . . , Xn sind unabhangig identisch verteilt, entweder gemaß ν0 oder gemaßν1, mit

H0 : Xi ∼ N (0, 1) =: ν0 ∀i H1 : Xi ∼ N (1, 1) =: ν1 ∀i.

Falls H0 zutrifft, mochten wir sie in hochstens 5% der Falle falschlicherweise ablehnen(mit einem Signifikanzniveau von 5%. Die Nullhypothese falschlicherweise abzulehnen(‘false positive’) bezeichnet man auch als Fehler erster Art). Damit durfen wir H0 nurablehnen, wenn

Z =X

σ/√n

=X

1/√n> 1.64.

Wir testen hier ‘einseitig’, denn nur große (positive) Werte von Z sprechen fur H1 undgegen H0. Es gilt:

WsH0(Z > 1.64) = 0.05.

Was passiert aber, wenn H1 zutrifft? Mit welcher Wahrscheinlichkeit lehnen wir dannH0 ab? Diese Wahrscheinlichkeit heißt die Macht eines Tests (testpower):

WsH1

(X

1/√n> 1.64

)= Ws

(X − 11/√n> 1.64−

√n

)= (1− Φ(1.64−

√n)),

wobei Φ die Verteilungsfunktion der Standardnormalverteilung ist.Fur n = 4 ist WsH1(|Z| > 1.64) ≈ 0.63, fur n = 10 etwa 0.87, fur n = 20 etwa 0.99.

18

Die Macht des z-Tests hangt ab von µ, σ, n und vom gewahlten Signifikanzniveau α(Ubung). Die Wahrscheinlichkeit, H0 falschlicherweise nicht abzulehnen (‘false negative’,1-testpower), bezeichnet man auch als Fehler zweiter Art. Insgesamt:

Entscheidung fur Entscheidung furH0 H1

mit W’keit mit W’keitH0 trifft zu 1− α α (Fehler 1. Art)H1 trifft zu β (Fehler 2. Art) 1− β

3 Der t-Test

Die z-Statistik hat den Vorteil, dass sie wegen des ZGWS fur alle Verteilungen von X(mit endlichem zweiten Moment) approximativ standardnormalverteilt ist. Aber der z-Test hat auch eine große Schwache: Er setzt voraus, dass σ bekannt ist. Dies ist aber imallgemeinen nicht der Fall; σ muss aus den Daten geschatzt werden.

3.1 Zur Geometrie von Strichprobenmittelwert und -varianz

3.1.1 Schatzer fur σ

Wegenσ2X = E

((X − µ)2

)bietet sich folgender Schatzer fur σ2

X an, falls µ bekannt ist

σ2X =

1n

∑i

(Xi − µ)2 (ist wg. E(X1 − µ)2 = σ2X erwartungstreu).

Aber µ ist selten bekannt und muss aus den Daten geschatzt werden: µ = X. Da aberX den Term

∑(Xi − c)2 minimiert, ist

1n

∑(Xi − X)2 ≤ 1

n

∑(Xi − µ)2.

Man kann zeigen, dass stattdessen

s2 =1

n− 1

n∑i=1

(Xi − X)2 (2)

ein erwartungstreuer Schatzer fur σ2 ist.

3.1.2 Geometrie des Mittelwertes

Wir fassen n Beobachtungen x1, . . . , xn ∈ R auf als einen Datenvektor X = (x1, . . . , xn) ∈Rn (Abb. 7). Sei D die ‘Hauptdiagonale’ in Rn, also

D := {c1|c ∈ R}, wobei 1 := (1, . . . , 1)T .

19

Abbildung 7: Darstellung von X = (x1, . . . , xn)T im Rn.

Dann ist x1 die orthogonale Projektion des Datenvektors X auf die Hauptdiagonale.Wir schreiben dafur

x1 = PDX.

Denn das Skalarprodukt

< 1,X− c1 >= 0 fur∑

xi − nc = 0, also: c = x.

Da durch die orthogonale Projektion c so bestimmt wird, dass die Lange der StreckeX − c1 minimal wird, bedeutet das gleichzeitig, dass

∑(xi − c)2 = ||X − c1||2 minimal

ist fur c = x. Anschaulich gesprochen zerlegen wir X in einen konstanten Anteil x1 undein Residuum R := X− x1:

X = x

1...1

+

x1 − x...

xn − x

= x1 +R,

wobei sowohl x1 als auch R als orthogonale Projektionen aufgefasst werden konnen:

x1 = PDX, und R = PD⊥X.

Abbildung 8: Geometrische Bedeutung von Mittelwert und Standardabweichung

20

3.1.3 Geometrie der geschatzten Standardabweichung

Will man den Vektor X ‘am besten’ (im Sinne des kleinsten euklidischen Abstandes)durch eine Zahl beschreiben, so wahlt man den Mittelwert x. Aber wie gut kann man X

durch x beschreiben? Wie weit ist X von der Hauptdiagonale D entfernt?

||R|| =√∑

(xi − x)2 =√n− 1 · s.

Je kurzer das Residuum R, desto naher ist X an D, desto weniger weichen die Beobach-tungen von ihrem Mittelwert ab.

3.2 Die t-Statistik

Ersetzt man in der z-Statistik σ durch s, so ist die resultierende Große (zumindest furkleine n) nicht mehr normalverteilt. Aber falls X1, . . . , Xn unabhangig und N (µ, σ2)-verteilt sind, hat sie eine bekannte Verteilung: die t-Verteilung mit n−1 Freiheitsgraden.

Satz: (William Gosset, alias ‘Student’, ∼ 1907)Seien X1, . . . , Xn unabhangig und N (µ, σ2)-verteilt. Dann hangt die Verteilung von

T :=X − µs/√n

nur von n ab (und z.B. nicht von σ). Genauer: T ist student-verteilt (t-verteilt) mit(n− 1) Freiheitsgraden (Beweis folgt).

3.2.1 Die t-Verteilung und ihre Verwandten

Die t-Verteilung ist glockenformig und (in aller Regel) um 0 zentriert. Sie gleicht derStandardnormalverteilung und geht fur n→∞ in sie uber. Fur kleine n hat sie jedochbreitere Schwanze (vgl. Figur 9, links): Eine t-verteilte ZVe T nimmt mit hoherer Wahr-scheinlichkeit extremere Werte an als eine normalverteilte ZVe Z. Die t-Statistik streutstarker als die z-Statistik, weil auch im Nenner eine Zufallsgroße steht. Das Gewichtaußerhalb des Intervalls [−1.96, 1.96] ist großer als 5%. Die neuen Schranken kann manin R mit der Funktion qt(0.025,df) berechnen (typische Werte siehe Tabelle 4).

df 2 5 10 15 30 ∞tkrit 4.30 2.57 2.23 2.13 2.04 1.96

Tabelle 4: Kritische Werte der t-Verteilung fur Ws(|T | > tkrit) = 0.05.

Die Dichte der t-Verteilung mit n Freiheitsgraden ist

ftn(x) =Γ((n+ 1)/2)√nπΓ(n/2)

(1 +

x2

n

)−n+12

, x ∈ R (vgl. Figur 9).

Eine t-verteilte ZVe kann beschrieben werden als Quotient von unabhangigen ZVen:

Sei T t(n)-verteilt, Z ∼ N (0, 1) und X ∼ χ2(n), so ist L(T ) = L(

Z√X/n

).

21

Normal− und t−VerteilungD

icht

e

−4 −2 0 2 4

0.0

0.1

0.2

0.3

0.4

df = 2df = 8

χχ2−Verteilung

Dic

hte

0 1 2 3 4 5 6 7

0.0

0.4

0.8 df = 1

df = 2df = 4df = 6

Abbildung 9: Die t- und die χ2-Verteilung

Die χ2-Verteilung Eine ZVe X ist χ2-verteilt mit n Freiheitsgraden, wenn

L(X) = L(Z2

1 + Z22 + · · ·+ Z2

n

), mit Z1, . . . , Zn unabhangig und N (0, 1)-verteilt.

Die Dichte der χ2-Verteilung mit n Freiheitsgraden ist

fχ2(n)(x) =1

Γ(n/2)12

(x2

)n2−1e−

x2 , x > 0 (vgl. Figur 9).

Mit EZ2i = 1 und Var(Z2

i ) = 2 sieht man: Ist X ∼ χ2(n), so ist E(X) = n undVar(X) = 2n. (Und am Rande sei bemerkt: An der Dichte sieht man sofort, dass dieχ2(2)-Verteilung der Exp(1/2)-Verteilung entspricht.)

Um zu zeigen, dass die T -Statistik t-verteilt ist, brauchen wir noch einige Voruber-legungen, die uns auch spater von Nutzen sein werden. Wir stellen dazu zunachst T aufandere Weise dar: Von X1, . . . , Xn unabhangig und N (µ, σ2)-verteilt gehen wir uber zu

Zi :=Xi − µσ

unabhangig und N (0, 1)-verteilt.

Dann ist

T =X − µs/√n

=X−µσ/√n

sσ

=√nZ

s/σ=

√nZ√

1n−1

∑(Zi − Z)2

.

Da Zi ∼ N (0, 1), ist auch√nZ ∼ N (0, 1). Es genugt also zu zeigen:

1.n∑i=1

(Zi − Z)2 ∼ χ2(n− 1)

2.√nZ und

n∑i=1

(Zi − Z)2 sind unabhangig.

22

3.3 Die Standardnormalverteilung auf Rn und der Satz von Fisher

Definition Z = (Z1, . . . , Zn)T heißt standardnormalverteilt auf Rn ⇐⇒ Z1, . . . , Znunabhangig und standardnormalverteilt in R. Damit ist die Wahrscheinlichkeitsdichtevon Z gegeben durch

ϕn(x) =n∏i=1

1√2πe−

x2i2 =

1(2π)n/2

e−Px2i2 =

1(2π)n/2

e−||x||2

2 , x = (x1, . . . , xn)T .

Insbesondere ist ϕn rotationssymmetrisch, weil es nur von der Lange ||x|| des Vektors x

abhangt.

Der Satz von Fisher Z sei standardnormalverteilt in Rn, H sei k-dimensionalerTeilraum von Rn. Dann gilt:

1. ||PHZ||2 ∼ χ2(k)2. PHZ und PH⊥Z sind unabhangig.

Zum Beweis brauchen wir ein

Lemma (uber die orthogonale Invarianz der Standardnormalverteilung)Z sei standardnormalverteilt in Rn. Dann gilt:

1. Fur M orthogonale n× n-Matrix (d.h. M ·MT = I) ist

Y := M · Z wieder standardnormalverteilt.

2. Zu jeder ONB von Rn sind die Koordinaten von Z wieder standardnormalverteiltund unabhangig.

Bew. 1.:Die Dichte von Y = MZ ist

g(y) =1

|detM |ϕn(M−1y), y ∈ Rn,

denn

Ws(M(z) ∈ B) = Ws(z ∈M−1(B)) =∫

M−1(B)

ϕn(z)dz

z=M−1(y)=

∫B

ϕn(M−1(y))1

|detM |dy.

Nun gilt aber

|detM | = 1 und ||M−1(y)||2 = yT (M−1)TM−1y = ||y||2,

und damit g(y) = ϕn(y).Bew. 2.:

23

Sei e1, . . . , en ONB von Rn und M := (e1, . . . , en)T . Dann sind die Koordinaten von Z

in der neuen ONB gegeben durch Y = MZ und damit nach 1. standardnormalverteiltund unabhangig. �

Beweis des Satzes von FisherWir wahlen eine ONB (e1, . . . , ek) vonH und erganzen sie zu einer ONBB := (e1, . . . , en)von Rn. Die Koordinaten von Z in B seien W1, . . . ,Wk,Wk+1, . . . ,Wn. Sie sind nach demvorausgegangenen Lemma wieder unabhangig und standardnormalverteilt. Mit

PHZ = W1e1 + . . .+Wkek = (W1, . . . ,Wk, 0, . . . , 0)T

PH⊥Z = Wk+1ek+1 + . . .+Wnen = (0, . . . , 0,Wk+1, . . . ,Wn)T

sind PHZ und PH⊥Z unabhangig, und ||PHZ||2 =∑k

i=1W2i ∼ χ2(k). �

Beweis, dass T ∼ t(n− 1) :Erinnerung: Wir hatten T folgendermaßen umgeschrieben:

T =X − µs/√n

=√nZ√

1n−1

∑(Zi − Z)2

,

mit Zi ∼ N (0, 1). Damit war√nZ ∼ N (0, 1), und es war zu zeigen:

1.n∑i=1

(Zi − Z)2 ∼ χ2(n− 1)

2.√nZ und

n∑i=1

(Zi − Z)2 sind unabhangig.

Da Z := Z1 Projektion von Z auf die Hauptdiagonale D ist und R := Z− Z Projektionvon Z auf das orthogonale Komplement D⊥, folgt die Behauptung. �

Man bemerke: Die Verteilung von T hangt damit nur von n ab (und nicht etwa von µund σ2)! Das ist ein Glucksfall, denn n ist der Stichprobenumfang und damit immerbekannt.

3.4 Anwendung: Eine Stichprobe

Wir wissen jetzt: Fur X1, . . . , Xn unabhangig und N (µ, σ2)-verteilt, mit unbekanntemµ und 0 < σ2 < ∞, ist die t-Statistik T =

√n(X − µ)/s nach t(n − 1)-verteilt. Was

nutzt diese Information in der Praxis?

1. Student-Konfidenzintervall Man kann (auch fur kleine n) ein exaktes Konfiden-zintervall I1−α(X) fur µ zum Niveau (1 − α) angeben. Wir schreiben fur das p-Quantilder t(n− 1)-Verteilung kurz qtn−1(p) := F−1

t(n−1)(p). Wegen

Wsµ [T ∈ [qtn−1(α/2), qtn−1(1− α/2)]] = 1− α

= Wsµ

[µ ∈

[X − s√

nqtn−1(1− α/2), X +

s√nqtn−1(1− α/2)

]]

24

ist

I1−α(X) :=[X − s√

nqtn−1(1− α/2), X +

s√nqtn−1(1− α/2)

]ein Konfidenzintervall fur µ zum Niveau (1− α).

2. Der Student-Test testet die Hypothese

Hθ : µ = θ.

Unter Hθ ist die Teststatistik

Tθ(X) :=√n(X − θ)

st(n− 1)-verteilt.

a) Klassischer t-Test mit Signifikanz α: Wegen

Wsθ [|Tθ(X)| ≥ qtn−1(1− α/2)] = α

lehne Hθ ab, falls |Tθ(X)| ≥ qtn−1(1 − α/2). Dann wird Hθ falschlicherweise abgelehntmit Wahrscheinlichkeit α.b) Allgemein gilt

Wsθ [|Tθ(X)| ≥ |τ |] = 2 · (1− Ft(n−1)(|τ |)).

Das heißt, nimmt |Tθ(X)| den Wert τ an, kann man Hθ ablehnen zum p-Wert 2 · (1 −Ft(n−1)(|τ |)).

Bemerkung Bei einem Test der Hypothese Hθ : µ = θ wird Hθ genau dann (zum p-Wert α) abgelehnt, wenn θ nicht in dem unter 1. hergeleiteten (1−α)-KonfidenzintervallI1−α(X) fur µ liegt.

3.5 Der zwei-Stichproben t-Test (fur gleiche Varianzen)

X1, . . . , Xm, Y1, . . . , Yn seien unabhangig, mit X1, . . . , Xm ∼ N (µ1, σ2), Y1, . . . , Yn ∼

N (µ2, σ2), und µ1, µ2, σ

2 seien unbekannt. Wie testet man die Hypothese

H0 : µ1 = µ2 (µ1 − µ2 = 0)?

Wir suchen eine geeignete Teststatistik (die nach Moglichkeit nicht von µ1, µ2, σ2 abhangt).

Bemerke:

X ∼ N (µ1, σ2/m), Y ∼ N (µ2, σ

2/n)

X − Y ∼ N(µ1 − µ2, σ

2

(1m

+1n

))N :=

X − Y − (µ1 − µ2)√1m + 1

n σ∼ N (0, 1)

Wir bauen daraus eine t-verteilte Statistik

T :=N√W/d

,

25

mit W ∼ χ2(d)-verteilt mit d = d(n,m) und N,W unabhangig. Wie geht das?Aus dem Satz von Fisher wissen wir:

1σ2

∑(Xi − X)2 ∼ χ2(m− 1) und unabhangig von X

1σ2

∑(Yi − Y )2 ∼ χ2(n− 1) und unabhangig von Y .

Die Vektoren X = (X1, . . . , Xm) und Y = (Y1, . . . , Yn) sind sowieso unabhangig, unddamit ist

W :=1σ2

[∑(Xi − X)2 +

∑(Yi − Y )2

]∼ χ2(n+m− 2)

und unabhangig von X, Y (und damit von N). Also gilt

T :=N√

W/(n+m− 2)∼ t(n+m− 2).

T lasst sich folgendermaßen darstellen:

T =X − Y − (µ1 − µ2)√

1m + 1

n

√1

n+m−2(∑

(Xi − X)2 +∑

(Yi − Y )2)=:

X − Y − (µ1 − µ2)√1m + 1

n sX,Y

,

mits2X,Y =

1n+m− 2

(∑(Xi − X)2 +

∑(Yi − Y )2

).

Konfidenzintervall fur µ1−µ2 Wie in der Situation einer Stichprobe ist ein (1−α)-Konfidenzintervall fur µ1 − µ2 gegeben durch die Grenzen

X − Y ±√

1m

+1nsX,Y qtm+n−2(1− α/2).

Test der Hypothese Hθ : µ1 − µ2 = θ Unter Hθ gilt

Tθ(X,Y) :=X − Y − θ√

1m + 1

n sX,Y

∼ t(m+ n− 2).

Damit konnen wir Hθ zum Niveau α ablehnen, falls

|Tθ(X,Y)| ≥ qtm+n−2(1− α/2).

Bemerkung Wieder wird bei einem Test die Hypothese Hθ genau dann (zum p-Wertα) abgelehnt, wenn θ nicht im (1− α)-Konfidenzintervall fur µ1 − µ2 liegt.

26

3.6 Zur Geometrie in der zwei-Stichproben-Situation

Wie in der Situation mit einer Stichprobe gibt es eine anschauliche geometrische In-terpretation der T -Statistik. Wir betrachten die Beobachtungen X1, . . . , Xm, Y1, . . . , Ynwiederum als Eintrage eines Datenvektors im Rm+n:

Z := (X1, . . . , Xm, Y1, . . . , Yn)t normalverteilt im Rm+n,

mit X1, . . . , Xm, Y1, . . . , Yn unabhangig und identisch verteilt gemaß N (µ, σ2). Wir be-trachten den 2-dimensionalen Unterraum M von Rm+n:

M := {(a, . . . , a︸︷︷︸m

, b, . . . , b︸︷︷︸n

)t|a, b,∈ R}.

M kann man in zwei orthogonale Unterraume zerlegen:

M := D ⊕ E,

wobei D wieder die Hauptdiagonale ist und E ihr orthogonales Komplement in M. Dund E werden jeweils aufgespannt durch die Einheitsvektoren

e0 :=1√

m+ n·

1...1

und e1 :=1√

1/m+ 1/n·

1/m...

1/m−1/n

...−1/n

.

Jetzt sehen wir, dass

||Z− PMZ||2 = ||PM⊥Z||2 =∑

(Xi − X)2 +∑

(Yj − Y )2 = (m+ n− 2) · s2X,Y

< Z, e1 > = ||PEZ|| = X − Y√1m + 1

n

und damit

T :=X − Y√

1m + 1

n · sX,Y

=||PEZ||

√m+ n− 2

||PM⊥Z||.

Die T -Statistik vergleicht die Langen von zwei orthogonalen Projektionen von Z auforthogonale Unterraume: E reprasentiert die Gruppenunterschiede zwischen X und Y

(die Variabilitat zwischen den Gruppen),M⊥ reprasentiert die Variabilitat innerhalb derGruppen. Wir werden dieselbe Zerlegung bei der klassischen Varianzanalyse antreffen,die auf analoge Art und Weise die Mittelwerte von mehr als zwei Gruppen miteinandervergleicht.

27

4 Tests fur gepaarte Stichproben

4.1 t-Test fur gepaarte Stichproben

Ein Beispiel aus Box-Hunter-Hunter, Statistics for Experimenters, Wiley, 1987:Zwei Materialien (A und B) fur Schuhsohlen sollen verglichen werden: Welches Materialnutzt sich weniger ab? Zehn Jungen erklaren sich als Testpersonen bereit, 2 Monate langnur die entsprechend besohlten Schuhe zu tragen.

• Versuchsplanung:

a) Wahle zufallig 5 der 10 Jungen und gebe ihnen Material A, die anderen be-kommen Material B. Vergleiche nach 2 Monaten die Abnutzung. Nachteil diesesVersuchsplans: Die Abnutzung der Materialien hangt nicht nur vom Material selbstab, sondern auch von den Individuen: Die Variabilitat kommt aus zwei verschie-denen Quellen, die entsprechenden Effekte sind konfundiert. Ein Materialeffekt istschwer zu entdecken, wenn die Unterschiede zwischen den Individuen groß sind.

b) Um Unterschiede zwischen den Individuen auszuschalten, geben wir jedem Jun-gen einen Schuh mit Sohle A und einen Schuh mit Sohle B (jeweils zufallig an denrechten oder linken Fuß).

Ergebnisse:

Tabelle 5: Bubenschuhe

Junge 1 2 3 4 5 6 7 8 9 10A 13.2 8.2 10.9 14.3 10.7 6.6 9.5 10.8 8.8 13.3B 14.0 8.8 11.2 14.2 11.8 6.4 9.8 11.3 9.3 13.6

• Auswertung:

a) Wir betrachten A1, . . . , A10, B1, . . . , B10 als unabhangig und normalverteilt undmachen einen t-Test fur zwei unabhangige Stichproben wie im letzten Abschnittbeschrieben. Das Ergebnis: p = 0.72. Wir konnen die Hypothese µA = µB alsonicht ablehnen - die beobachteten Unterschiede konnen in mehr als 70% aller Ver-suche durch Zufall zustande kommen.

b) Wir haben in a) schon wieder nicht berucksichtigt, dass die Daten eine Paa-rungsstruktur haben, dass also jeder Wert Ai genau einem Wert Bi entspricht.Die Variabilitat zwischen den Individuen ging in a) noch in die T-Statistik ein.Wir konnen aber leicht davon Gebrauch machen, dass sich bei manchen Jungendie Sohlen grundsatzlich weniger abnutzen als bei anderen: die Stichproben sindpositiv korreliert (siehe Abb. 10). Wir betrachten die Differenzen:

Di := Ai −Bi

und untersuchen mit einem t-Test fur eine Stichprobe die Nullhypothese, dass dieVerteilung dieser Differenzen Mittelwert Null hat. Das Ergebnis: p = 0.008.

28

●

●

●

●

●

●

●

●

●

●

6 8 10 12 14

68

1012

14A

B

Abbildung 10: Gepaarte Stichproben

4.1.1 Voraussetzungen und Spezialfall

Die einzige Voraussetzung beim t-Test fur gepaarte Stichproben ist, dass der Effekt (indem Fall die Unterschiede in der Abnutzung des Materials) normalverteilt ist und un-abhangig fur jedes Individuum. Man kann also den Effekt E einer Behandlung auf einMerkmal untersuchen, ohne die Verteilung des Merkmals in der untersuchten Popula-tion uberhaupt zu kennen. Allein der Effekt, der sich in der Differenz Di = Ai − Biaußert, muss normalverteilt sein, nicht die zugrundeliegenden Verteilungen von A undB. Dennoch sind es interessante Spezialfalle, wenn A und B normalverteilt sind:

Normalverteilung Im Fall A1, . . . , An, B1, . . . , Bn unabhangig und normalverteilt,mit gleicher Varianz und EA = µA, EB = µB wurde man fur die Hypothese µA = µBeinen Test fur zwei unabhangige Stichproben durchfuhren, mit

Tungepaart :=A− B√

2n sA,B

,

wobei T t-verteilt ware mit 2n− 2 Freiheitsgraden.Was passiert, wenn man statt dessen eine Paarungsstruktur einbaut und die DifferenzenDi = Ai − Bi mit einem Ein-Stichproben-Test betrachtet? (Das kann man machen,denn die Differenzen von unabhangigen Normalverteilten sind ja wieder normalverteilt).Der wesentliche Unterschied besteht darin, dass man Freiheitsgrade einbußt, denn dieentsprechende Teststatistik

Tgepaart :=A− B√1n sA−B

ist unter der Nullhypothese t-verteilt mit nur (n − 1) Freiheitsgraden. Damit verliertman ein wenig Testmacht, was fur große n allerdings vernachlassigbar ist. Ansonstenunterscheiden sich die beiden Teststatistiken im Prinzip nicht, denn fur unabhangigeAi, Bi mit gleicher Varianz σ2 ist Var(A−B) = σ2

A−B = 2σ2A.

Sind dagegen Ai und Bi positiv korreliert (wie im Beispiel der Bubenschuhe), so istVar(A − B) = 2σ2

A − 2Cov(A,B). Infolgedessen wird der Nenner von Tgepaart kleiner,die Teststatistik wird großer, und man gewinnt an Testmacht.

29

Positiv korrelierte Normalverteilte Modell:A1, . . . , An seien unabhangig undN (µ, σ2)-verteilt, also Ai = µ+σZi fur Zi ∼ N (0, 1). B1, . . . , Bn seien unabhangig und normalver-teilt. Aber Ai und Bi seien nicht unabhangig, sondern positiv korreliert in der folgendenForm:

Bi = Ai + δ + αZ ′i, mit δ ∈ R, α > 0, Z ′i ∼ N (0, 1).

Dann gilt:

Bi −Ai = δ + αZ ′i, und damit

Var(Di) = α2

4.2 Wilcoxon-Test fur gepaarte Stichproben/ Vorzeichenrangtest

Im t-Test war die einzige benotigte Modellannahme die Normalverteilung der Differen-zen:

L(D) = N (δ, σ2), bzw. L(D − δ) = N (0, σ2).

Der Wilcoxon-Test kommt ohne Normalverteilung aus. Er testet die Hypothese Hµ, dassdie Verteilung der Differenzen symmetrisch um µ ist. (Alternativen sind hier also sowohl:Eine symmetrische, aber nicht um µ zentrierte Verteilung, oder eine Verteilung, die nichtsymmetrisch ist.) Die Teststatistik verwendet die mit den Vorzeichen versehenen Range:

Wµ(D) =n∑i=1

1{Di>µ}Rang(|Di − µ|).

Unter Hµ ist Wµ so verteilt wie

S =n∑i=1

i · Vi,

wobei Vi unabhangige, faire {0, 1}−Munzwurfe sind. Die Verteilung von S kennt R unterdem Namen signrank (rsignrank, psignrank, etc.). Die Verteilung von S ist symme-trisch um

ES =∑

iEVi =12

∑i =

n(n+ 1)4

.

Ein Test auf Hµ : δ = µ liefert als p-Wert entsprechend die Wahrscheinlichkeit, dass dieTeststatistik S um mehr als |Wµ − n(n+1)

4 | von n(n+1)4 abweicht. In R einfach mit dem

Befehl: wilcox.test(A-B), oder wilcox.test(A,B,paired=TRUE). Im obigen Beispiel:p = 0.014

4.3 Vorzeichentest

Eine weitere simple Statistik verwendet ausschließlich die Vorzeichen der Differenzen:Unter der Nullhypothese dass Ai und Bi dieselbe Verteilung haben, erwarten wir ge-nauso viele positive wie negative Differenzen Di = Ai − Bi. Die Anzahl der positivenDifferenzen:

S2 := |{i|Di > 0}|ist somit binomialverteilt mit Parametern n und 1/2. Im obigen Beispiel: Nur 2 von 10Differenzen sind positiv. Der p-Wert: 2 pbinom(2,10,0.5)≈ 0.11. Der Vorzeichentestist weniger machtig als der Vorzeichenrangtest, weil er nicht berucksichtigt, dass diepositiven Differenzen vom Betrag her kleiner sind als die negativen Differenzen.

30

4.4 Permutationstest von Fisher

Und noch eine Statistik fur gepaarte Stichproben:

S3 :=n∑i=1

(sgn Di) ·Di.

Unter der Nullhypothese gilt: Gegeben die Werte der Differenzbetrage |Di|, ist S3 soverteilt wie

∑|Di|Vi, wobei Vi wieder unabhangige faire {−1, 1}-Munzwurfe sind.

5 Vergleich von mehr als 2 Stichproben und das normalelineare Modell

5.1 Vom t-Test zur einfaktoriellen Varianzanalyse

5.1.1 Geometrie beim t-Test

Wir erinnern uns an das Vorgehen beim t-Test zum Vergleich von zwei unabhangigenStichproben und an die geometrische Situation in Sektion 3.6. Gegeben sind zwei Stich-proben X11, . . . , X1n1 und X21, . . . , X2n2 , alle Xij unabhangig und normalverteilt mitgleicher Varianz σ2, aber eventuell mit unterschiedlichen Erwartungswerten E[X1i] = µ1

fur i = 1, . . . , n1, und E[X2j ] = µ2 fur j = 1, . . . , n2. Das Modell ist also:

Xij = µi + σZij , wobei Zij ∼ N (0, 1) und unabhangig, i = 1, 2, j = 1, . . . , ni,

mit µ1, µ2 ∈ R. Oder in Vektorschreibweise:

X = µ+ σZ,

mit X = (X11, . . . , X1n1 , X21, . . . , X2n2)t, µ = (µ1, . . . , µ1, µ2, . . . , µ2)t und Z standard-normalverteilt im Rn1+n2 . Das allgemeine Modell sagt also

µ ∈M mit M := {(µ1, . . . , µ1, µ2, . . . , µ2)t|µ1, µ2 ∈ R}, dim(M) = 2.

Wir untersuchen die Nullhypothese H0 : µ1 = µ2 und damit die Frage, ob µ eigentlichin einem eindimensionalen Unterraum D von M liegt mit

D := {(µ, . . . , µ)t|µ ∈ R}.

Um dies zu tun, zerlegen wir den Vektor X in die orthogonalen Projektionen auf dieUnterraume D, E (mit M = D ⊕ E) und auf M⊥:

X = PDX + PEX + PM⊥X

= x

1...1...

+

x1 − x

...x2 − x

...

+

x11 − x1

...x21 − x2

...

,

wobei xi := 1ni

∑nij=1Xij den Gruppenmittelwert in Gruppe i bezeichnet und x :=

1n1+n2

∑i,j Xij den globalen Mittelwert von X.

31

Geometrische Interpretation der t-StatistikNun ist die t-Statistik gerade:

T :=x1 − x2

s√

1/n1 + 1/n2

=||PEX||

||PM⊥X||/√n1 + n2 − 2

.

Unter der Nullhypothese: µ ∈ D ist wegen D⊥E: PEµ = 0, und damit PEX = PEσZ.Damit ist nach dem Satz von Fisher unter der Nullhypothese die t-Statistik (saloppgeschrieben):

T =

√χ2(1)√

χ2(n1 + n2 − 2)/√n1 + n2 − 2

=N (0, 1)√

χ2(n1 + n2 − 2)/√n1 + n2 − 2

∼ t(n1+n2−2).

Definition (Fisher-Verteilung)Seien X1 ∼ χ2(k1) und X2 ∼ χ2(k2) unabhangige Zufallsgroßen. Dann heißt die Vertei-lung von

X1/k1

X2/k2

Fisher-Verteilung mit k1 und k2 Freiheitsgraden (man sagt auch k1 Zahler- und k2

Nennerfreiheitsgrade).

BeobachtungIst T student-verteilt mit n Freiheitsgraden, so ist T 2 Fisher -verteilt mit 1 und n Frei-heitsgraden. Das Quadrat der obigen t-Statistik ist also Fisher(1, n1+n2−2)-verteilt. Ein(zweiseitiger) t-Test zum Test der Hypothese µ1 = µ2 gegen µ1 6= µ2 zum Niveau 1− αmit Hilfe von entsprechenden student(n1 +n2−2)-Quantilen ist also aquivalent zu einemTest, der die Statistik T 2 und die entsprechenden Quantile der Fisher(1, n1 + n2 − 2)-Verteilung verwendet.

5.1.2 Verallgemeinerung auf k(> 2) Gruppen

Wie kann man untersuchen, ob die Mittelwerte von mehr als 2 Gruppen sich unter-scheiden? Ein Beispiel: 24 Wunden von Ratten wurden mit je einer von 4 verschiedenenBehandlungen (A, B, C, D) behandelt. Hat die Behandlung einen Einfluss auf die Blut-gerinnungszeit? Abbildung 11 zeigt, dass die Mittelwerte der vier Gruppen sich unter-scheiden. Hatte das Zufall sein konnen? Wie groß sind die Unterschiede zwischen denvier Mittelwerten, die wir durch Zufall erwarten?

Wir gehen ganz analog zum 2-Stichproben-t-Test vor. Das zugrundeliegende Modellfur k Gruppen:

X1j = µ1 + σZ1j fur j = 1, . . . , n1

...Xkj = µk + σZkj fur j = 1, . . . , nk,

alsoX = µ+ σZ mit µ ∈M

32

Behandlung

Ger

innu

ngsz

eit (

s)

A B C D

6065

70

Abbildung 11: Blutgerinnungszeit bei vier verschiedenen Wundbehandlungen

undM := {(µ1, . . . , µ1, . . . , µk, . . . , µk)t|µ1, . . . , µk ∈ R}, dim(M) = k. Wir untersuchenwieder die Nullhypothese µ1 = . . . = µk, also

µ ∈ D mit D := {(µ, . . . , µ)t|µ ∈ R}

und projizieren dazu X orthogonal auf die zueinander orthogonalen Unterraume D, E(mit M = D ⊕ E und somit dim(E) = k − 1) und M⊥:


= x

1...1...

+

x1 − x

...xk − x

...

+

x11 − x1

...xk1 − xk

...

.

Wir betrachten die Teststatistik

F :=||PEX||2/(k − 1)||PM⊥X||2/(n− k)

. (3)

In jedem Fall gilt µ ∈ M, und damit ist PM⊥µ = 0, und es gilt nach dem Satz vonFisher

1σ2||PM⊥X||2 = ||PM⊥Z||2 ∼ χ2(n− k).

Fur den Zahler von F gilt unter der Nullhypothese (µ ∈ D) weiterhin PEµ = 0 unddamit

1σ2||PEX||2 = ||PEZ||2 ∼ χ2(k − 1).

Also gilt unter der Nullhypothese, d.h. wenn alle Gruppenmittelwerte gleich sind:

F :=||PEX||2/(k − 1)||PM⊥X||2/(n− k)

∼ Fisher(k − 1, n− k).

Trifft die Nullhypothese aber nicht zu, d.h. ist µ /∈ D, so ist PEµ 6= 0 und damit derZahler von F tendenziell großer.

33

Klassische einfaktorielle VarianzanalyseGroße Werte von F deuten also auf Abweichungen von der Nullhypothese hin. Ein Testder Nullhypothese zum Niveau 1− α hatte damit den Annahmebereich

{F |F ≤ qf(1− α, k − 1, n− k)},

wobei qf(1−α, k−1, n−k) das 1−α-Quantil der Fisher-Verteilung mit k−1 und n−kFreiheitsgraden bezeichnet.

Die klassische Varianzanalyse (kurz: ANOVA - ANalysis Of VAriance) zerlegt alsoeigentlich nur die Variabilitat der Daten, ||X||2 in einen gemeinsamen Anteil, ||PDX||2,in einen Anteil, der die Unterschiede zwischen den Gruppenmittelwerten reprasentiert,||PEX||2, und in einen Anteil, der die zufallige Variabilitat innerhalb der Gruppen re-prasentiert, ||PM⊥X||2. Wegen E||PM⊥X||2 = σ2(n− k) ist der Nenner von F ,

s2 := ||PM⊥X||2/(n− k),

ein erwartungstreuer Schatzer fur σ2. Wegen Pythagoras ist

||X||2 = ||PDX||2 + ||PEX||2︸︷︷︸Treatmentquadratsumme

+ ||PM⊥X||2︸︷︷︸Fehlerquadratsumme

= nx2 +∑i

ni(xi − x)2 +k∑i=1

ni∑j=1

(xij − xi)2.

Im Blutgerinnungsbeispiel haben wir n1 = 4, n2 = n3 = 6, n4 = 8, x1 = 61, x2 =66, x3 = 68, x4 = 61, n =

∑ni = 24 und x = 64 und damit folgende Zerlegung:


mit xij = x+ (xi − x) + (xij − xi),

626063596367...

59

=

646464646464

...64

+

−3−3−3−3

22...−3

+

1−1

2−2−3

1...−2

,

mit den entsprechenden Quadratsummen ||PEX||2 = 228 und ||PM⊥X||2 = 112, unddamit

F =||PEX||2/(k − 1)||PM⊥X||2/(n− k)

=228/3112/20

≈ 13.6.

Abbildung 12 zeigt, dass so große Unterschiede unter der Nullhypothese sehr seltenvorkommen (p < 0.001).

34

0 5 10 15

0.0

0.2

0.4

0.6

x

Dic

hte

13.6

Abbildung 12: Dichte der Fisher-Verteilung mit 3 und 20 Freiheitsgraden.

5.2 Das normale lineare Modell

Die Varianzanalyse ist ein Beispiel einer statistischen Fragestellung im Rahmen eineslinearen Modells. Allgemein haben wir es zu tun mit einem Modell der Form

X = µ+ σZ,

wobei Z standardnormalverteilt auf Rn, σ ≥ 0 und µ ∈ M, mit M k-dimensionalerTeilraum des Rn.

Im Rahmen dieses Modells konnen wir µ und σ schatzen und lineare Hypothesentesten:

Schatzung von µ und σ2

Fur den kleinste-Quadrate-Schatzer µ von µ gilt:∑

(xi − µi)2 = minµ∈M∑

(xi − µi)2,und daher

||X− µ||2 = minm∈M

||X−m||2,

und so ist der kleinste-Quadrate-Schatzer wieder die orthogonale Projektion von X aufM:

µ = PMX = PMµ+ PMσZ = µ+ PMσZ.

Also ist µ erwartungstreuer Schatzer fur µ. Analog mit R := PM⊥X = σPM⊥Z:

s2 :=||R||2

n− k= σ2 ||PM⊥Z||2

n− k.

Wegen ||PM⊥Z||2 ∼ χ2(n − k) ist s2 erwartungstreuer Schatzer fur σ2, und nach demSatz von Fisher unabhangig von µ.

µ ist auch der ML-Schatzer fur µ, denn fur X = (X1, . . . , Xn) unabhangig undXi ∼ N (µi, σ2) ist die Likelihoodfunktion gegeben durch∏

i

1√2πσ2

exp−(xi − µi)2σ2

=1

(2πσ2)n/2exp−

∑(xi − µi)2

2σ2.

Sie wird maximal, wenn die Summe der Residuenquadrate∑

(xi − µi)2 minimal wird,also fur den Kleinste-Quadrate-Schatzer µ.

35

5.2.1 Testen linearer Hypothesen

Wir gehen ganz analog zur Varianzanalyse vor, nur wahlen wir als ModellraumM einenbeliebigen k-dimensionalen Teilraum von Rn und als Hypotheseraum D einen beliebigend-dimensionalen Teilraum vonM. Wir untersuchen wie gehabt die Nullhypothese µ ∈ D.Wie bisher zerlegen wir den Datenvektor X durch orthogonale Projektion auf D undM(mit M = D ⊕ E):

X = PDX + PEX + PM⊥X.

Wie bisher ist unter der Nullhypothese µ ∈ D und damit

PEX = σPEZ =⇒ 1/σ2||PEX||2 ∼ χ2(dim(E)).

Analog giltPM⊥X = σPM⊥Z =⇒ 1/σ2||PM⊥X||2 ∼ χ2(n− k),

und damit ist unter der Nullhypothese die Teststatistik

F :=||PEX||2/(k − d)||PM⊥X||2/(n− k)

Fisher(k − d, n− k)-verteilt.

Abbildung 13: Normales Lineares Modell in der Nussschale

5.3 Paarweise Vergleiche, multiples Testen

In der klassischen einfaktoriellen Varianzanalyse wird die Nullhypothese untersucht, dassalle Gruppenmittelwerte gleich sind. Welche Schlusse erlaubt ein signifikantes Ergebnis?Welche Gruppen unterscheiden sich?

Zunachst kann man zuruckkehren zu den paarweisen Vergleichen, die wir aus dem t-Test kennen. Fur ein fest gewahltes Paar (j, `) von Gruppen konnen wir obiges Vorgehen

36

einschranken, um nur Unterschiede zwischen diesen beiden Gruppen zu untersuchen: Wirwahlen von M einen zweidimensionalen Teilraum

Mj,` := {(0, . . . , 0, a, . . . , a, 0, . . . , 0, b, . . . , b, 0, . . . , 0)t|a, b ∈ R},

bei dem die Eintrage aller anderen Gruppen Null sind. Analog

Dj,` := {(0, . . . , 0, a, . . . , a, 0, . . . , 0, a, . . . , a, 0, . . . , 0)t|a ∈ R} ⊂ Mj,`,

und Ej,` sei das orthogonale Komplement von Dj,` in Mj,`. Wie bisher ist

PEj,`X = (0, . . . , 0, xj , . . . , xj , 0, . . . , 0, xl, . . . , xl, 0, . . . , 0)t,

und damit gilt unter der Nullhypothese µj = µ`:

Fj,` :=||PEj,`X||2/1

||PM⊥X||2/(n− k)∼ Fisher(1, n− k),

alsoTj,` :=

√Fj,` ∼ student(n− k).

Also ist

Ij,` :=

[xj − xl −

√1nj

+1n`s · qt(1− α/2), xj − xl +

√1nj

+1n`s · qt(1− α/2)

]ein Konfidenzintervall fur µj − µ` zum Niveau 1 − α. Man bemerke: Dieses Vorgehenist eine Mischung aus klassischem t-Test und klassischer einfaktorieller Varianzanalyse,da zwar nur zwei Gruppen verglichen werden, aber der Nenner der Teststatistik, s2, ausdem gesamten Datensatz geschatzt wird.

5.3.1 Bonferroni-Korrektur

Hat man ein bestimmtes Paar (j, `) im Auge, eignet sich das eben beschriebene Vorgehengut. Aber was, wenn man alle Paare miteinander vergleichen will? Testet man jedes der(k2

)Paare jeweils zum Signifikanzniveau α, so ist naturlich die Wahrscheinlichkeit, dass

mindestens ein Test anschlagt, obwohl µ1 = µ2 = . . . = µk, nicht mehr α, sonderngroßer. Bei

(k2

)unabhangigen Tests ware sie 1 − (1 − α)(

k2), und fur 0 < α < 1 ist

(1− α)(k2) < 1− α.

Eine einfache, aber recht konservative Moglichkeit, dieses Problem in den Griff zubekommen, ist, α durch α/

(k2

)zu ersetzen und die entsprechenden Quantile fur den Test

zu verwenden. Die Idee dahinter:

Ws(mind. ein Test schlagt an) = Ws({Test 1 schlagt an} ∪ . . . ∪ {Test(k

2

)schlagt an})

≤(k2)∑i=1

Ws(Test i schlagt an) =(k

2

)α

Also ist 1−(k2

)α eine einfache obere Schranke fur (1− α)(

k2). Allerdings sind die resul-

tierenden Konfidenzintervalle recht groß.

37

5.3.2 Simultane Konfidenzintervalle nach Tukey

Wenn alle Gruppen gleich groß sind (n1 = n2 = . . . = nk =: m), kann man direktsimultane Konfidenzintervalle fur alle Paare µj − µ` bilden. Sei dafur zunachst

Icj,` :=[xj − x` − c · s/

√m, xj + x` + c · s/

√m]

∀j, `

Fur jedes Paar (j, `) bilden wir so ein Intervall mit derselben Breite, die von c abhangt.Wie groß ist c zu wahlen, damit

1− α = Ws(µj − µ` ∈ Icj,` ∀(j, `))= Ws(xj − x` − c · s/

√m ≤ µj − µ` ≤ xj − x` + c · s/

√m ∀(j, `))

= Ws(max(j,`)

√m/s |(xj − µj)− (x` − µ`)| ≤ c)

= Ws(σ

s

(maxj

√m(xj − µj)

σ−min

`

√m(x` − µ`)

σ

)≤ c)

Die Terme√m(xi−µi)

σ , i = 1, . . . , k, sind unabhangig und N (0, 1)-verteilt. Das bedeutet,wir konnen fur k unabhangige und N (0, 1)-verteilte ZVe Z1, . . . , Zk die obige Wahr-scheinlichkeit umschreiben zu

1− α = Ws

(maxj Zj −min` Z`√

s2/σ2≤ c

).

Weiterhin gilt

L(s2

σ2

)= L

(||PM⊥X||2

(n− k)σ2

)=

1n− k

χ2(n− k),

und s und√m(xi−µi)

σ sind unabhangig. Daher muss man c so wahlen, dass

1− α = Ws

(maxj Zj −min` Z`√

Y/g≤ c

),

wobei Z1, . . . , Zk unabhangig und N (0, 1)-verteilt und unabhangig von Y ∼ χ2(g). DieVerteilung von

maxj Zj −min` Z`√Y/g

heißt studentsche Spannweitenverteilung mit Parametern k (Anzahl der Gruppen) undn− k (Anzahl der Freiheitsgrade im Nenner).

ZahlenbeispieleFur k = 4, n = 24, α = 0.05 findet man c =qtukey(1 − α, k, n − k) = 3.96. Mit derBonferroni-Korrektur hatte man (mit

√1/nj + 1/n` =

√2/m) ein Konfidenzintervall

mit Faktor qt(1− α/12, 6 + 6− 2) ·√

2 ≈ 3.28 ·√

2 = 4.63.Fur k = 6, n = 36, α = 0.05 erhalt man mit Tukey: c = 4.3, mit Bonferroni 5.41.Fur k = 10, n = 60, α = 0.05 erhalt man mit Tukey: c = 4.68, mit Bonferroni 6.39.

38

5.4 Lineare Kontraste

Außer paarweisen Tests kann man auch andere Fragen stellen, z.B.: Ist die Differenzµ1 − µ2 so groß wie die Differenz µ3 − µ4, also: µ1 − µ2 − (µ3 − µ4) = 0? All dies sindBeispiele fur sogenannte lineare Kontraste, d.h. Ausdrucke der Form

k∑i=1

ciµi mit∑

ci = 0.

Wir suchen nach simultanen Konfidenzintervallen fur alle linearen Kontraste von µ,d.h. fur alle

κc =∑

ciµi mit∑

ci = 0.

Konfidenzintervall fur einen KontrastFur einen einzigen Kontrast κc =

∑ciµi gilt mit κc :=

∑cixi und τ2

c :=∑c2i :

L(κ) = L(∑

cixi

)= N

(∑ciµi, σ

2/m · τ2c

).

Also ist ein (1− α)–Konfidenzintervall fur κ gegeben durch

Iκ :=[κ− s/

√m · τc · qt(1− α/2, n− k), κ+ s/

√m · τc · qt(1− α/2, n− k)

].

(Wir schreiben κ statt κc.)

Simultane Konfidenzintervalle fur alle linearen KontrasteWir suchen ein γ so, dass fur alle c = (c1, . . . , ck) mit

∑ci = 0 gilt:

Iγκ :=[κ− s/

√m · τc · γ, κ+ s/

√m · τc · γ

]ist ein (1− α)-Konfidenzintervall fur den Kontrast κ, d.h. so, dass

1− α = Ws(κ ∈ Iγκ ∀κ) = Ws{

maxc

||κ− κ||τc

√m

s≤ γ

}(4)

Dazu machen wir zunachst zwei Beobachtungen:

Beobachtung 1Fur µ ∈M = D ⊕ E zerlegen wir µ in seine orthogonalen Projektionen auf D und E:

µ = PDµ+ PEµ =: δ + β.

E wird dadurch beschrieben, dass alle x ∈ E senkrecht auf (1, . . . , 1)t stehen, also:

E := {(a1, . . . , ak)t|∑

ai = 0}.

Wegen∑ci = 0 ist daher c ∈ E und damit

κc =∑

ciµi =< c, µ >

=< c,PDµ > + < c,PEµ >=< c,PEµ >=∑

ciβi,

mit β = (β1, . . . , βk)t. Außerdem definiert jedes Element von E einen Kontrast.

39

Beobachtung 2Wegen κc =

∑ciβi ist

||κ− κ|| = |k∑i=1

ci(βi − βi)|C.−S.≤

√∑c2i︸︷︷︸

=τc

√∑(βi − βi)2.

Bei Cauchy-Schwarz (C.-S., < x, y >≤ ||x||·||y||) gilt genau dann Gleichheit, wenn x undy kollinear sind, d.h. oben ist κ− κ maximal, wenn es ein a so gibt, dass ci = a(βi− βi)∀i = 1, . . . , k, d.h. wenn die ci proportional zu (βi−βi) sind. Ein solches a gibt es, dennwegen β ∈ E und β ∈ E ist auch (β − β) ∈ E. Wir kehren zuruck zu Gleichung (4) undsehen:

maxc

(||κ− κ||

τc

)=√∑

(βi − βi)2,

und damit 1− α = Ws{

maxc

||κ− κ||τc

√m

s≤ γ

}= Ws(m/s2 ·

∑(βi − βi)2 ≤ γ2)

= Ws(

σ2||PEZ||2

σ2||PM⊥Z||2/(n− k)≤ γ2

)= Ws

(||PEZ||2/(k − 1)||PM⊥Z||2/(n− k)

≤ γ2

k − 1

).

Fur ein simultanes Konfidenzintervall zum Niveau (1 − α) fur alle Kontraste κc mussalso (mit n = m · k) gelten:

γ :=γ2

k − 1= qf(1− α, k − 1, k(m− 1)),

also γ ist das 1− α-Quantil der Fisher(k − 1, k(m− 1))-Verteilung, und damit sind mit

γ =√

(k − 1)γ =√

(k − 1)qf(1− α, k − 1, k(m− 1))

die IntervalleIγκ :=

[κ− s/

√m · τc · γ, κ+ s/

√m · τc · γ

]simultane (1− α)-Konfidenzintervalle fur alle Kontraste κ.

5.5 Kruskal-Wallis-(H)-Test

Wie auch beim t-Test gibt es auch fur den Wilcoxon-Test ein analoges nichtparametri-sches Verfahren zum Vergleich von k Stichproben, den Kruskal-Wallis-Test. Analog zumWilcoxon-Test werden nicht die Zahlen, sondern die Range verwendet.

40

Erinnerung: Wilcoxon-TestGegeben seien zwei unabhangige Stichproben

X11, . . . , X1n1 unabhangig und verteilt gemaß π1


Der Wilcoxon-Test testet die Nullhypothese, dass beide Stichproben aus derselben Ver-teilung stammen, also dass π1 = π2, gegen die Alternative, dass eine Verteilung gegendie andere verschoben ist. Die Teststatistik war die Summe der Range in Gruppe 1:

S :=n1∑j=1

R1j ,

mit ES = n1(n1 + n2 + 1)/2, Var(S) = n1n2(n1 + n2 + 1)/2, und fur n1, n2 →∞

S − ES√Var(S)

→ N (0, 1) in Verteilung.

Verallgemeinerung auf k Stichproben Das Modell und die Hypothesen sind ana-log:


......

Xk1, . . . , Xknk unabhangig und verteilt gemaß πk

Sei n :=∑ni. Wir untersuchen die Nullhypothese:

H0 : π1 = . . . = πk.

Wieder ersetzen wir dazu die Zahlen Xij durch ihre Range Rij aus den Gesamtdaten.Der mittlere Rang aller n Zahlen ist

R =1n

n(n+ 1)2

=n+ 1

2.

Wir untersuchen die Abweichungen der mittleren Range in den k Gruppen, Ri :=1ni

∑nij=1Rij , i = 1, . . . , k vom mittleren Rang R :

H :=12

n(n+ 1)

k∑i=1

ni(Ri − R)2.

Sofern die πi kontinulierlich sind, sind unter der Nullhypothese die Range R11, . . . , Rknkeine zufallige Permutation von 1, . . . , n, und es gilt fur n→∞ und ni/n→ pi > 0

H −→ χ2(k − 1) in Verteilung.

Diese Approximation ist (als grobe Faustregel) fur k = 3, ni ≥ 5 bzw. k ≥ 4 und ni ≥ 4akzeptabel.

41

BemerkungHaufig schreibt man H in anderer Form:

H =12

n(n+ 1)

∑ S2i

ni− 3(n+ 1),

wobei Si :=∑

j Rij die Summe der Range in Gruppe j bezeichnet.

6 Lineare Regression

6.1 Modell und Parameterschatzer

Bei der linearen Regression wird der lineare Zusammenhang zwischen mehreren Varia-blen untersucht. Die klassische Fragestellung: Gegeben seien je n Beobachtungen zweierVariablen X und Y :

X = (x1, . . . , xn), Y = (y1, . . . , yn).

Wenn X und Y in annahernd linearer Beziehung stehen, vermuten wir eine Beziehungder Form

yi = β0 + β1x+ Fehler.

Wir bestimmen zunachst die Koeffizienten β0 und β1:

Kleinste QuadrateIn der klassischen Fragestellung werden die Beobachtungen xi und yi als fest betrachtet,und wir bestimmen β0 und β1 nach der Methode der kleinsten Residuenquadrate:

n∑i=1

(yi − (β0 + β1xi))2 = min!

Die Gerade {(x, y)|y = β0 + β1x} heißt Regressionsgerade.

Maximum LikelihoodLineare Regression kann man aber auch im Rahmen des normalen linearen Modellsbetrachten, und damit sind wieder die ML-Schatzer identisch mit den Kleinste-Quadrate-Schatzern. Diesmal werden nur die xi als feste reelle Zahlen betrachtet:

Yi = β0 + β1xi + σZi,

wobei Z1, . . . , Zn unabhangig und N (0, 1)-verteilt, also Yi ∼ N (β0 + β1xi, σ2), oder

Y = β01 + β1X + σZ.

Das heißt, die Likelihoodfunktion ist gegeben durch

L =n∏i=1

1√2πσ2

exp−(Yi−(β0−β1xi))

2

2σ2 =1

(2πσ2)n/2exp−

∑(Yi − (β0 + β1xi))2

2σ2.

Wiederum wird L maximiert fur minimale Residuenquadratsumme∑

(Yi−(β0+β1xi))2,das heißt, die ML-Schatzer und die Kleinste-Quadrate-Schatzer fur β0 und β1 stimmenuberein.

42

Alternative Schreibweise: Das lineare RegressionsmodellWir betrachten das Vorhergehende im Lichte des normalen linearen Modells:

Y = β01 + β1X + σZ = µ+ σZ,

mitµ ∈M := {α01 + α1X|α0, α1 ∈ R}, dim(M) = 2.

Alternativ:

Y =

1 x1...

...1 xn

︸︷︷︸

=:C

(β0

β1

)+ σZ.

Die n×2-Matrix C nennt man Systemmatrix, β = (β0, β1)t die Regressionskoeffizienten.Die Modellannahme ist µ ∈M := {Cβ|β ∈ R2}.Die kleinste-Quadrate-Schatzer fur β in Matrixschreibweise:

(Y− Cβ)⊥Cβ =⇒ (Y− Cβ)tC = 0

⇐⇒ YtC = βtCtC

⇐⇒ CtY = CtCβ

Dies ist eindeutig losbar genau dann, wenn CtC regular ist, also vollen Rang hat, d.h.wenn die Spaltenvektoren von C linear unabhangig sind.

Die Schatzer fur die univariate RegressionIn der klassischen Fragestellung Y = β0 + β1X+ Fehler haben die Parameterschatzerfolgende anschaulichen Eigenschaften.

1) Der Schwerpunkt (x, y) der Datenwolke liegt auf der Regressionsgeraden:

(Y− Cβ)⊥1 =⇒∑

(yi − β0 − β1xi) = ny − nβ0 − β1nx = 0.

2) Die Steigung β1 hangt ab von σx, σy und der Korrelation zwischen x und y:(Y− Cβ)⊥X =⇒

0 =∑

(yi − β0 − β1xi)xi (und wegen∑

(yi − β0 − β1xi)x = 0 gilt:)

=∑

(yi − β0 − β1xi)(xi − x)

=∑

(yi − β1xi)(xi − x)

=∑

((yi − y)− β1(xi − x))(xi − x)

=⇒ β1 =∑

(yi − y)(xi − x)∑(xi − x)2

=Cov(x, y)

s2x

= rx,ysysx.

Falls r = 1, ist β1 = sysx

: wachst x um sx, so wachst y um sy. Falls |r| < 1, so fallt dieSteigung entsprechend geringer aus. Falls r = 0, so ist die Steigung der Regressionsge-raden immer 0, denn die Kenntnis von x liefert keine Informationen uber die Große vony.

43

6.2 Die Geometrische Bedeutung von r und der Multiple Korrelati-onskoeffizient

Um den Korrelationskoeffizienten fur die multiple Regression verallgemeinern zu konnen,machen wir zunachst die folgende Beobachtung: Es gilt

r2 =||PEY||2

||PD⊥Y||2(5)

Bedeutung: Welcher Anteil der Variabilitat in Y (||PD⊥Y||2) wird durch die Vorhersagedurch X (||PEY||2) eingefangen (’erklart’)? Beweis von (5): Setze

M := {β01 + β1X|β0, β1 ∈ R}D := {β01|β0 ∈ R}E := {c(X− x1)|c ∈ R}

Wir beobachten:1) PEY = c(X− x1) mit c =

P(yi−y)(xi−x)P

(xi−x)2. Warum?

< Y− c(X− x1),X− x1 >= 0

=⇒∑

yi(xi − x)− c∑

(xi − x)2 = 0

=⇒ c =∑yi(xi − x)∑(xi − x)2

=∑

(yi − y)(xi − x)∑(xi − x)2

,

denn∑y(xi − x) =

∑yxi − nyx = 0. Damit ist

||PEY||2 = c2∑

(xi − x)2 =(∑

(yi − y)2(xi − x)2)2∑(xi − x)2

.

2) Wegen PD⊥Y = Y− y1 folgt

||PD⊥Y||2 =∑

(yi − y)2.

Also gilt:

r2 =||PEY||2

||PD⊥Y||2=

∑((yi − y)(xi − x))2∑

(xi − x)2∑

(yi − y)2

Der multiple KorrelationskoeffizientWir betrachten nun ein Modell, in dem Y nicht nur durch eine Variable X vorhergesagtwerden soll, sondern durch viele Variablen X1, . . . , X`−1:

Y = Cβ + σZ = β0 + β1X1 + β2X2 + . . .+ β`−1X`−1 + σZ= µ+ σZ

mit µ ∈M := {Cβ|β ∈ R`}. Wir stellen die Frage: Helfen die Variablen X1, . . . , X`−1 beider Vorhersage von Y? Oder ist die Verbesserung der Vorhersage gegenuber dem Modell

44

Y = β0 + σZ nicht großer als durch Zufall zu erwarten, d.h. wenn PEY = PEσZ? Umdiese Frage zu untersuchen, konnte man

R2 :=||PEY||2

||PD⊥Y||2

verwenden. Unter der Nullhypothese gilt fur den Zahler

1/σ2||PEY||2 = ||PEZ||2 ∼ χ2(`− 1)

und fur den Nenner

1/σ2||PD⊥Y||2 = ||PD⊥Z||2 ∼ χ2(n− 1).

Da aber sowohl Zahler als auch Nenner die Projektion auf E enthalten und nicht mehrunabhangig sind, verwendet man

||PEY||2

||PM⊥Y||2=

||PEY||2||PD⊥Y||2

||PM⊥Y||2||PD⊥Y||2

=R2

1−R2,

denn ||PM⊥Y||2 = ||PD⊥Y||2 − ||PEY||2. Und hier sieht man leicht die Verwandtschaftzur bekannten F -verteilten Teststatistik: Unter der Nullhypothese gilt

R2/(`− 1)(1−R2)/(n− `)

∼ F (`− 1, n− `).

Fur ` = 2 folgt damit sofort:

r2 · (n− 2)1− r2

∼ F (1, n− 2)

bzw.r√n− 2√

1− r2∼ t(n− 2). (6)

Der Term in (6) ist die klassische Statistik zum Test der Korrelation normalverteilterZufallsvariablen.

6.2.1 Ein Konfidenzintervall fur die Korrelation

Wir betrachten hier nur den Fall von normalverteilten Zufallsvariablen. Zunachst unter-scheiden wir zwischen der ’wahren’ Korrelation ρ zweier Zufallsvariablen X und Y undder empirischen Korrelation r ihrer Realisierungen (x1, . . . , xn), (y1, . . . , yn). Anhandvon r mochten wir ein Konfidenzintervall fur ρ herstellen. Dies ist allerdings nicht leicht,da die Verteilung und Varianz von r stark von ρ abhangen (Abb. 14). Wahrend fur ρ ≈ 0die Verteilung, zumindest fur großere n, annahernd normal ist, ist sie fur |ρ| >> 0 (undnicht genugend große n) deutlich schief.

45

ρρ = 0, n=100

r

Dic

hte

−1.0 −0.5 0.0 0.5 1.0

01

23

4

ρρ = 0, n=10

r

Dic

hte

−1.0 −0.5 0.0 0.5 1.0

0.0

0.4

0.8

1.2

ρρ = 0.5, n=100

r

Dic

hte

−1.0 −0.5 0.0 0.5 1.0

01

23

45

ρρ = 0.5, n=10

r

Dic

hte

−1.0 −0.5 0.0 0.5 1.0

0.0

0.5

1.0

1.5

Abbildung 14: Verteilungen von r fur ρ = 0 (links) und ρ = 0.5 (rechts), n = 100 (oben)und n = 10 (unten).

−1.0 −0.5 0.0 0.5 1.0

−4

−2

02

4

r

z

Abbildung 15: Fisher’s z-Transformation

46

Ein Ausweg: Fisher’s z-Transformation

z :=12

log(

1 + r

1− r

).

Diese Große streckt Werte von r nahe 1 (Abb. 15) und ist fur normalverteilte Zufalls-variablen annahernd N (0, 1/

√n− 3)-verteilt. Die Vorteile: 1) Die Standardabweichung

hangt nicht mehr von ρ ab, sondern nur von n, und muss daher nicht geschatzt werden.2) Die Form der Verteilung hangt auch nicht mehr von ρ ab.

Ein Konfidenzintervall fur ρ baut man dann einfach dadurch, dass man ein KIfur den z-transformierten Wert erstellt und die Grenzen rucktransformiert. Beispiel:n = 10, r = 0.9.

=⇒ z = 1.472. 95%−KI : z ± 1.96 ∗ 1/√

7 ≈ z ± 0.741 : [0.831, 2.213].

Rucktransformation ergibt ein 95%-Konfidenzintervall fur ρ von [0.624,0.976]. Offen-sichtlich ist dieses nicht symmetrisch um r, aber entspricht der Schiefe der Verteilungund liegt im Bereich der moglichen Werte fur ρ.

Ein Vergleich zwischen r1 und r2 kann analog durchgefuhrt werden, in dem r1 undr2 z-transformiert werden, dann ist approximativ

z1 ∼ N (0, 1/(n1 − 3)) und z2 ∼ N (0, 1/(n2 − 3)),

und damit gilt approximativ

z1 − z2√1

n1−3 + 1n2−3

∼ N (0, 1).

6.3 Ein Konfidenzbereich fur β

Definition Kovarianzmatrix Sei X := (X1, . . . , Xk)t eine Rk-wertige Zufallsvariablemit Erwartungswert µ = EX. Dann heisst die Matrix mit Eintragen

cij := Cov(Xi, Xj) = E[(Xi − µi)(Xj − µj)]

die Kovarianzmatrix von X.

Definition Normalverteilung auf Rk X heißt normalverteilt auf Rk mit Mittelwertµ und Kovarianzmatrix C (schreibe: X ∼ N (µ,C)) genau dann wenn

EX = µ, Cov(Xi, Xj) = cij∀i, j

und ∀a ∈ Rk ist aTX normalverteilt auf R1.

47

Behauptung:Sei wie im normalen linearen Regressionsmodell Y = Cβ + σZ, mit Z standardnormal-verteilt auf Rn. Dann ist β normalverteilt mit Eβ = β und Kovarianzmatrix σ2(CTC)−1,denn:

β = (CTC)−1CTY = (CTC)−1CT (Cβ + σZ)

= (CTC)−1CTCβ + (CTC)−1CTσZ

= β + (CTC)−1CTσZ

Wegen EZ = 0 ist β erwartungstreu, und

Var((CTC)−1CTσZ) = (CTC)−1CTσ2C(CTC)−1 = σ2(CTC)−1.

Satz (von Cochran) (z.B. Shao Mathematical Statistics 2nd edition, 2003, S. 27)Sei V ∼ N (0,Σ) normalverteilt auf Rk und habe Σ vollen Rang k. Dann ist

V TΣ−1V ∼ χ2(k).

Was bedeutet dies fur den Konfidenzbereich fur β? Mit β − β ∼ N (0, σ2(CTC)−1)bedeutet dies

(β − β)T · CTC

σ2· (β − β) ∼ χ2(`),

also ist ein Konfidenzbereich fur β zum Niveau 1− α gegeben durch

β + {v ∈ R`|vT · CTC

σ2· v ≤ qchisq(1− α, `)}.

6.4 Multivariate Regression

Wir haben bereits oben Modell und Parameterschatzer fur das normale lineare Regres-sionsmodell mit mehr Variablen besprochen: Y = β0 + β1X1,+β2X2 + . . . + β`−1X`−1,mit β = (CTC)−1CTY und gezeigt, dass mit Σ := σ2(CTC)−1

β ∼ N (β,Σ).

6.4.1 Orthogonale Regressoren

Im Fall von orthogonalen Regressoren Xi⊥Xj ∀i 6= j kann man M in orthogonalelineare Unterraume zerlegen, die von den Xi aufgespannt werden. Die βi erhalt mandurch orthogonale Projektion von Y auf diese Unterraume. Diese Projektionen sindunabhangig, und damit sind die Schatzer βi unabhangig, d.h. Σ ist eine Diagonalmatrix,und es gilt

βi =1

||Xi||2< Y, Xi >= βi +

σ

||Xi||< Z,

Xi

||Xi||> =⇒ βi ∼ N (βi,

σ2

||Xi||2).

Aus der Unabhangigkeit folgt auch, dass ∀i die Schatzer von βi dieselben sind, egalwie viele der Variablen X1, . . . , Xi−1, Xi, . . . , X`−1 zur Vorhersage von Y herangezogenwerden. (Bei nicht orthogonalen Regressoren ist das idR nicht der Fall)

48

Beispiel Lineare Regression mit X = 0.

7 Schatzverfahren

In der Regel gibt es fur einen Parameter mehrere Schatzer. Welchen soll man nehmen?Welche Schatzverfahren gibt es, und wie bewertet man die Schatzer?

Beispiel: X1, . . . , Xn seien unabhangig und uniformverteilt auf [0, `]. Als Schatzer fur` bieten sich an:

ˆ1 := 2X, denn EX = `/2.

ˆ2 :=

n+ 1n

maxi=1,...,n

(Xi).

Beide Schatzer sind erwartungstreu: Im zweiten Fall kompensiert der Faktor n+1n dafur,

dass max(Xi) =: M < `. Um zu zeigen, dass Eˆ2 = `, berechnen wir zunachst die

Verteilungsfunktion von M :

F (x) = Ws(M ≤ x) =(x`

)n=⇒ d

dxF (x) =

n

`nxn−1, 0 ≤ x ≤ `.

=⇒ EM =

`∫0

x( n`nxn−1

)dx =

n

n+ 1`.

Welchen der beiden Schatzer sollte man nun nehmen? Wir vergleichen die erwartetenquadratischen Fehler:

E(

(ˆ1 − `)2

)EX=`= Var(2X) =

4n

Var(Xi)

=4n

(E(X2

i )− E(Xi)2)

=4n

`∫0

x2

`dx− `2

4

=4n· `

2

12

=1

3n· `2.

E(

(ˆ2 − `)2

)Eˆ

2=`= Var(n+ 1n

M

)=

(n+ 1)2

n2E(M2)− `2.

Mit

E(M2) =

`∫0

x2( n`2xn−1

)dx =

n

n+ 2· `2

folgt

E(

(ˆ2 − `)2

)= `2

[(n+ 1)2

n(n+ 2)− 1]

=1

n(n+ 2)· `2.

Der mittlere quadratische Fehler von ˆ2 fallt mit n−1, aber der mittlere quadratische

Fehler von ˆ2 fallt sogar mit n−2.

49

7.1 Der Mittlere quadratische Fehler

Warum der mittlere quadratische Fehler? Der mittlere quadratische Fehler isteiner der gangigsten Maßstabe zur Bewertung von Schatzern. In der Tat liefert er beimVergleich von Schatzern dieselben Resultate wie jede ‘vernunftige’ vergleichbare Gute-funktion EL: Hat von zwei Schatzern θ1, θ2 eines Parameters θ1 den kleineren quadrati-schen Fehler, so ist auch EL(θ1) < EL(θ2).Fur eine ‘vernunftige vergleichbare’ Gutefunktion L gelte

• L(θ) = 0: Schatzt man den wahren Parameter, so ist der Fehler Null.

• L ist zweimal differenzierbar und hat an der Stelle θ ein lokales Minimum, so dassL′(θ) = 0 und L′′(θ) > 0.

Dann kann man L in der Nahe von θ beschreiben durch c · (θ − θ)2:

L(θ) = L(θ)︸︷︷︸=0

+L′(θ)︸︷︷︸=0

·(θ − θ) +L′′(θ)

2︸︷︷︸=c(θ)>0

·(θ − θ)2 + . . .

Ist θ nahe bei θ, kann man die weiteren Terme vernachlassigen. Damit

EL(θ) ≈ c(θ)E((θ − θ)2).

Hat also θ1 einen kleineren mittleren quadratischen Fehler als θ2:

E((θ1 − θ)2)� E((θ2 − θ)2), so ist auch EL(θ1)� EL(θ2).

Konsistenz und der mittlere quadratische Fehler Grob gesprochen, bezeichnetdie Konsistenz eines Schatzers die Eigenschaft, dass die Abweichung des Schatzers vomzu schatzenden Wert mit wachsendem Stichprobenumfang gegen Null geht. Formal:Ein Schatzer θn heißt konsistenter Schatzer fur θ, falls fur alle ε > 0 gilt:

Ws(|θn − θ| > ε)→ 0 fur n→∞.

Fur die Konsistenz eines Schatzers reicht es, wenn sein mittlerer quadratischer Fehlergegen Null geht, also wenn

E(

(θn − θ)2)→ 0 fur n→∞.

Um dies einzusehen, stellen wir fest: |θn−θ|2 ≥ ε2 gilt auf der Menge A := {|θn−θ| > ε}.Damit konnen wir schreiben

|θn − θ|2 ≥ ε2 · IA,

wobei IA = 1 die Indikatorvariable der Menge A ist. Wir bilden die Erwartungswerte:

E(

(θn − θ)2)≥ ε2 ·Ws(|θn − θ| > ε).

Fur festes ε und n→∞ geht die linke Seite gegen Null und damit auch die rechte. �

50

Zerlegung des mittleren quadratischen Fehlers Falls θ erwartungstreu ist, soentspricht der mittlere quadratische Fehler gerade der Varianz von θ:

E(

(θ − θ)2)

Eθ=θ= Var(θ).

Ist θ aber nicht erwartungstreu, d.h.

b(θ) = E(θ)− θ 6= 0 (b...‘bias’),

so kann man den mittleren quadratischen Fehler zerlegen in die Varianz und das Quadratder Verzerrung:

E(

(θ − θ)2)

= E(

(θ − Eθ + Eθ − θ)2)

= E(

(θ − Eθ)2)

+ 2b(θ)E(

(θ − Eθ))

+ b2(θ)

= Var(θ) + b2(θ).

Haufig sind Var(θ) und b(θ) gegenlaufig - je kleiner das eine, desto großer das andere.Es kann daher manchmal sinnvoll sein, eine kleine Verzerrung in Kauf zu nehmen, umden erwarteten quadratischen Fehler zu minimieren:

Kontrollierte Verzerrung eines Schatzers. Ein Beispiel Wir wollen die Wahr-scheinlichkeit p fur das Eintreten eines Ereignisses A schatzen und fuhren dazu n un-abhangige Versuche durch. Seien I1, . . . , In mit Ii = 1, falls A eingetreten ist, und Nullsonst.Ein plausibler Schatzer fur p ist die relative Haufigkeit der Erfolge:

p1 :=1n

∑i

Ii.

p1 ist erwartungstreu und hat mittleren quadratischen Fehler E[(p1 − p)2

]= Var(p1) =

1nVar(I1) = p(1−p)

n . Aber: Dieser Schatzer hat nicht zwangslaufig minimalen quadrati-schen Fehler!Betrachte einen zweiten Schatzer

p2 := γp1 fur geeignet zu wahlendes γ.

Mit E(p2) = γp ist p2 fur γ 6= 1 nicht erwartungstreu. Aber:

E[(p2 − p)2

]= Var(p2) + (γp− p)2 =

1n

[γ2p(1− p) + np2(1− γ)2

].

Abhangig von p und γ kann dies kleiner sein als p(1−p)n . Fur γ = 0.9 ist z.B.

E[(p1 − p0)2

]= E

[(p2 − p0)2

]fur p0 ≈

11 + 0.05n

.

Fur p < p0 hat p2 kleineren quadratischen Fehler, fur p > p0 ist p1 zu bevorzugen. Wennalso bekannt ist, dass p sicher kleiner ist als ein festes p0, so kann man γ geeignet wahlen(Ubung).

51

Aus dem letzten Beispiel haben wir gelernt, dass die Gute eines Schatzers auch vomWert des Parameters abhangen kann: Fur einen Parameterbereich kann ein Schatzerbesser sein als ein anderer, fur einen anderen Parameterbereich kann es anders aussehen.Es gibt in der Regel keine Schatzmethode, fur die der mittlere quadratische Fehler furalle Werte von θ kleiner ware als fur alle anderen Methoden. Aber es gibt eine Methode,die meistens gut funktioniert:

7.2 Das Maximum-Likelihood Prinzip

SeienX1, . . . , Xn u.i.v. gemaß einer Verteilung, die von θ abhangt. Gegeben (X1, . . . , Xn) =(x1, . . . , xn), ist der Maximum-Likelihood-Schatzer von θ diejenige Zahl θML, die demEreignis {(X1, . . . , Xn) = (x1, . . . , xn)} maximales Gewicht zuordnet.Im diskreten Fall ist θML so definiert, dass WsθML

((X1, . . . , Xn) = (x1, . . . , xn)) maxi-mal wird.Hat X1 Dichte fθ(x), so soll

∏ni=1 fθ(xi) maximal werden.

Bsp. 1: Munzwurf mit Erfolgsparameter p. Von n Versuchen beobachten wir k Er-folge, mit

Ws(X = k) =(n

k

)pk(1− p)n−k.

Das Maximum dieser Funktion findet man am einfachsten durch Ableiten ihres Loga-rithmus, denn da log(x) eine streng monotone Funktion ist, ist das Maximum von log fan derselben Stelle wie das Maximum von f .

log(Ws(X = k)) = c+ k log p+ (n− k) log(1− p) und damitddp

log(. . .) =k

p− n− k

1− p.

(Fur p = 0 oder p = 1 gilt Ws(X = 0) = 1 bzw. Ws(X = n) = 1. In beiden Fallen stimmtder ML-Schatzer pML mit p uberein.) Fur die anderen Falle findet man durch Nullsetzender Ableitung den ML-Schatzer pML = k/n, die relative Haufigkeit der Erfolge.

Bsp. 2: Uniformverteilung Seien X1, . . . , Xn unabhangig und Unif[0, `]-verteilt. Ge-sucht ist der ML-Schatzer fur `. Die Dichte von X ist

f`(x) =

{1/` x ∈ [0, `]0 sonst

=⇒ f`(x1, . . . , xn) =

{1/`n xi ∈ [0, `]∀ i0 sonst.

Figur 16 zeigt den Graphen von f`(x1, . . . , xn): Er ist maximal bei ` = max(xi), dochfindet man das Maximum hier nicht durch Ableiten, weil f` dort nicht differenzierbarist.

7.2.1 Eigenschaften von Maximum-Likelihood-Schatzern

1.: Ist θ ML-Schatzer von θ, so ist der ML-Schatzer einer (injektiven) Funktion h(θ)gerade h(θ). Warum? Wenn das Ereignis {X = x} fur θML maximales Gewicht hat, sohat es auch unter h(θ) = h(θML) maximales Gewicht.

52

f l((x 1

,, ……,, x

n))maxxi

1/`n

`

Abbildung 16: f`(x1, . . . , xn) fur X1, . . . , Xn unabhangig und Unif[0, `]-verteilt.

Beispiel Die Lebensdauer X eine Sorte Gluhbirnen sei Exp(λ)-verteilt, d.h.

Ws(X ≤ x) = 1− e−λx.

Wir suchen die Wahrscheinlichkeit, dass eine Gluhbirne weniger als 200 Stunden brennt,und beobachten dazu die Lebensdauern x1, . . . , xn von n Gluhbirnen.

• Der ML-Schatzer von λ istλML =

1X,

wie man leicht durch Logarithmieren und Ableiten der gemeinsamen Dichte fλ =∏i λe−λxi = λne−λ

Pxi feststellt.

• Der ML-Schatzer von Ws(X ≤ 200) = 1− e−λ200 ist damit 1− e−λML200.

• Da die Exponentialverteilung nur von λ abhangt, ist auch 1−e−λMLx ML-Schatzerder zugehorigen Verteilungsfunktion, und fur jedes Intervall [a, b] ist der ML-Schatzer von Ws(X ∈ [a, b]) gegeben durch e−λa − e−λb.

• Es kommt nicht auf die Parametrisierung an: Bei (im englischen Sprachraum ubli-cher) Parametrisierung der Exponentialverteilung mit β = 1/λ:

λe−λx → 1βe− xβ gilt βML =

1

λML

= X.

2.: Fur eine (im folgenden naher beschriebene) bestimmte Klasse von Problemen istdie ML-Methode asymptotisch (d.h. fur n → ∞) im Sinne des mittleren quadratischenFehlers mindestens so gut wie jede andere Schatzmethode. (ML ist asymptotisch effizi-ent)

Die Problemklasse ‘glatter’ Probleme Sei R der Wertebereich der ZVe X, d.h.

R :=

{{x|Ws(X = x) > 0} fur X diskrete ZVe,kleinstes Intervall [a, b], mit

∫ ba f(x)dx = 1 fur X stetige ZVe.

Es soll gelten

53

• fθ(x) (bzw. Wsθ(X = x)) hat stetige zweite Ableitung in θ,

• Die Ableitungen mussen auf sinnvolle Weise beschrankt sein: Fur jedes kleineIntervall J von θ-Werten soll∫

maxθ∈J

∣∣∣∣ ∂∂θfθ(x)∣∣∣∣dx <∞ und

∫maxθ∈J

∣∣∣∣ ∂2

∂2θfθ(x)

∣∣∣∣ dx <∞• Fur alle θ soll maxx∈R |fθ(x)| <∞.

Dies trifft fur die meisten gangigen Verteilungen zu. Ausnahmen bilden die uniformeVerteilung und die doppelte Exponentialverteilung.

Die Fisher-Information Hat X Dichte fθ(x), die von θ abhangt, so ist die Fisher-Information der Verteilung von X in θ definiert als

I(θ) = E

[(∂

∂θlog fθ(x)

)2].

Also:

I(θ) =

{∑i

(∂∂θ log Wsθ(X = xi)

)2Wsθ(X = xi) im diskreten Fall,∫ (

∂∂θ log fθ(x)

)2fθ(x)dx im Dichtefall.

Beispiel X ∼ Exp(λ). Dann ist I(θ) = 1/λ2, denn

log fλ(x) = log λ− λx, ∂

∂λlog fλ(x) =

1λ− x

=⇒ I(λ) =∫ (

1λ− x)2

fλ(x)dx =1λ2.

In glatten Problemen gibt es eine untere Schranke fur den mittleren quadratischenFehler:

Satz Sei ϑn ein konsistenter Schatzer fur ϑ in einem Problem der oben beschriebenenKlasse. Dann gilt:

E[(ϑn − ϑ)2] ≥ 1nI(ϑ)

+ o

(1n

)(Beweisidee in Breiman, Statistics: With a view toward applications, S. 94)

Satz Der ML-Schatzer ϑMLn sei eindeutig (d.h. die ML-Gleichungen haben eine ein-

deutige Losung). Dann gilt

E[(ϑMLn − ϑ)2] =

1nI(ϑ)

+ o

(1n

)(ohne Beweis)

54

3.: ML-Schatzer sind in der Regel konsistent (bis auf pathologische Falle).

4.: Asymptotische Normalitat des ML-Schatzers Unter (relativ schwachen) Regula-ritatsvoraussetzungen gilt die schwache Konvergenz:

L(√n(ϑML

n − ϑ))→ N (0, I(ϑ)−1) (n→∞).

Die Regularitatsvoraussetzungen sind haufig erfullt: Lϑ(x) := log fϑ(x) muss zweimalstetig differenzierbar sein, und es muss gelten∫

∂2

∂ϑ2fϑ(x)dx = 0. (7)

Was impliziert das anschaulich? Unter dieser Bedingung entspricht die Fisher-Informationgerade dem Erwartungswert der zweiten Ableitung von L.Wir schreiben kurz: f ′ fur ∂

∂ϑf und f ′′ fur ∂2

∂ϑ2 f .

I(ϑ) = E[L′ϑ(X)2] = E[f ′ϑ(X)2

f2ϑ(X)

]= −E

[fϑ(X)f ′′ϑ(X)− f ′ϑ(X)2

f2ϑ(X)

]= −E[L′′ϑ(X)],

denn mit∫

∂2

∂ϑ2 fϑ(x)dx = 0 ist E(f ′′ϑ(X)/fϑ(X)) = 0. Die Fisher-Information entsprichtalso der (negativen) mittleren Krummung: Je starker Lϑ(x) im Mittel an der Stelle ϑgekrummt ist, desto kleiner ist die asymtotische Varianz des ML-Schatzers.

Man rechnet Bedingung (7) leicht nach. Z.B. fur X exponentialverteilt und entspre-chend fϑ(x) = ϑe−ϑx:

∂2

∂ϑ2fϑ(x) =

∂

∂ϑ

(e−ϑx − xϑe−ϑx

)= −xe−ϑx − x

(e−ϑx − xϑe−ϑx

)= −2xe−ϑx + x2ϑe−ϑx∫

∂2

∂ϑ2fϑ(x)dx = − 2

ϑEX + Var(X) + (EX)2 = − 2

ϑ2+

1ϑ2

+1ϑ2

= 0

7.3 Die Momentenmethode

Eine weitere gangige Methode, Schatzer zu bestimmen, ist die Methode der Momente:Man vergleiche die k-ten Momente der Verteilung:

EXk =∫xkfϑ(x)dx

mit den Momenten der Stichprobe:

1n

∑xki ,

k = 1, 2, . . . ,. Haufig stimmen Momentenschatzer und ML-Schatzer uberein:Munzwurf mit Erfolgsparameter p: E( 1

n

∑Ii) = p =⇒ p = I .

55

Exponentialverteilung X1, . . . , Xn unabh., Exp(λ)-verteilt: Mit EX = 1λ folgt λ = 1

X.

Nicht so bei der Uniformen Verteilung auf [0, `], wie wir bereits gesehen haben: DerMomentenschatzer von ` ist ˆ = 2X wegen EX = `

2 . Der ML-Schatzer ist M =max(X1, . . . , Xn).Momentenschatzer sind nach Konstruktion unverzerrt, und konsistent, und sie sindmanchmal leichter zu ermitteln als ML-Schatzer. Aber sie haben eventuell einen ho-hen mittleren quadratischen Fehler.

7.4 Kleinste Quadrate

Wir mochten einen Datensatz moglichst gut mit einer bestimmten Funktion V (ϑ) be-schreiben. Einfachstes Beispiel: X1, . . . , Xn sollen ’moglichst gut’ mit einer konstantenFunktion V (ϑ) = ϑ1 beschrieben werden, und zwar so, dass die Quadratsumme derResiduen zwischen V und X minimal wird:∑

(ϑ−Xi)2 = min .

Den Schatzer in diesem Fall kennen wir schon: ϑ = X. Da minimale Quadratsummendurch orthogonale Projektionen entstehen, haben wir auch im t-Test solche Schatzerkennen gelernt. Wir werden weitere Beispiele in der Varianzanalyse und der linearenRegression sehen.

Haufig stimmen die Schatzer vieler Methoden uberein: Fur X1, . . . , Xn unabhangigund identisch normalverteilt mit Mittelwert µ ist X der kleinste-Quadrate-Schatzer, derMomentenschatzer und der ML-Schatzer fur µ.

7.5 Bootstrap-Konfidenzintervalle

Ein Schatzer hangt von den zufalligen Beobachtungen ab und ist daher variabel. Oftkann man seine exakte Varianz oder gar seine exakte Verteilung bestimmen, wie bei dert-Statistik. Wenn das nicht moglich ist, gibt es evtl. asymptotische Aussagen (z.B. dieapproximative Normalitat der t-Statistik fur unabhangige und identisch verteilte ZVen).

Was kann man tun, wenn man die Varianz eines Schatzers weder exakt noch asym-ptotisch ermitteln kann? Man konnte sie empirisch bestimmen, indem man den Versuchoft wiederholt. Z.B.: Ein Versuch mit n unabhangigen p-Munzwurfen Ii. Wir schatzenp = 1/n

∑Ii (und wissen: Var(p) = p(1− p)/n). Ohne dieses Wissen konnten wir durch

haufiges Wiederholen dieses Versuches auf empirische Weise die Varianz schatzen.Aber das Wiederholen eines Versuches ist oft teuer oder unmoglich. Daher legt die

Bootstrapmethode die empirische Verteilung statt der wahren Verteilung zu Grunde,und es werden ‘Bootstrapstichproben’ aus der empirischen Verteilung gezogen: n Werte,unabhangig (d.h. mit Zurucklegen) und identisch verteilt (gemaß der empirischen Ver-teilung). Die Idee ist, dass die Varianz des Schatzers in den Bootstrapstichproben dieVarianz des Schatzers aus den echten Stichproben approximiert.

Wir haben also Beobachtungen X = (X1, . . . , Xn), unabhangig und identisch verteiltgemaß einer Verteilung ν. Wir schatzen ein Merkmal m(ν) der Verteilung als Funktionder Daten: h(X). Wie ist die Verteilung von h(X)−m(ν)? Statt neue Stichproben aus derOriginalverteilung ν zu ziehen, ziehen wir Stichproben X∗ der Große n (mit Zurucklegen)

56

aus der empirischen Verteilung νX. Idee:

L(h(X)−m(ν)) ≈ L(h(X∗)−m(νX))

Die Verteilung von h(X∗)−m(νX) wird in der Praxis durch Monte-Carlo-Simulationermittelt, die obere und untere Grenze des Konfidenzintervalls fur m(ν) schatzt manaus den empirischen Quantilen der durch Simulation ermittelten Verteilung.

Ubrigens ist das Vorgehen im Munzwurffall ganz ahnlich, obwohl wir die Formel furdie Varianz kennen: Wir schatzen Var(p) = p(1 − p)/n mit p(1 − p)/n und ermittelndamit die Varianz des Schatzers fur den Fall p = p. Das ist dasselbe als wurde manBootstrapstichproben aus der empirischen Verteilung mit Parameter p ziehen.

7.6 Suffizienz und die Minimierung d. mittleren quadratischen Fehlers

Wir beginnen mit einem Beispiel:X1, . . . , Xn seien unabhangige {0, 1}−Munzwurfe mit Erfolgsparameter p ∈ (0, 1). K :=∑Xi sei die Anzahl der Erfolge. Dann hangt die Verteilung von (X1, . . . , Xn) fur festes

K = k (d.h. gegeben das Ereignis {K = k}) nicht mehr von p ab:

Wsp({(X1, . . . , Xn) = (x1, . . . , xn)}{|∑

Xi = k})

=Wsp({(X1, . . . , Xn) = (x1, . . . , xn)} ∩ {K = k})

Wsp({K = k})

=pk(1− p)n−k(nk

)pk(1− p)n−k

=1(nk

) .Anschaulich: Uber p konnen wir nicht mehr Informationen gewinnen als durch Betrach-tung der Funktion K =

∑Xi der Daten. Die genauen Zeitpunkte, zu denen die K

Erfolge auftreten, hangen von p nicht mehr ab. Die Statistik K nennt man deswegenauch suffiziente Statistik fur p:

Definition Eine Statistik V (X) heißt suffizient fur ϑ, falls die bedingte Verteilung vonX, gegeben V (X), nicht von ϑ abhangt.

Beispiel: Multinomialverteilung Y1, . . . , Yn seien unabhangig und nehmen Wertein k Kategorien an, mit Ws(Yi in Kategorie j) = pj ∀i, j = 1, . . . , k,

∑pj = 1. Der

interessierende Parameter ist also der Vektor der Gewichte ~p = (p1, . . . , pk). Dann istder Vektor der Besetzungszahlen X = (X1, . . . , Xk) suffizient fur ~p (Rechnung analogzum obigen Beispiel).

Jede Dichte kann man mit Hilfe einer suffizienten Statistik V (X) geeignet faktorisie-ren in eine Funktion von V (X), die von ϑ abhangt, und einen Faktor, der nicht von ϑabhangt:

Der Faktorisierungssatz von Neyman und Pearson Die Verteilung von X besitzedie Dichte (bzw. Gewichte) fϑ(x). Dann sind aquivalent:

57

i) V (X) ist suffiziente Statistik fur ϑ.ii) Die Dichte zerfallt in eine Faktorisierung der Gestalt

fϑ(x) = gϑ(V (x)) · h(x),

wobei die Funktion h nicht von ϑ abhangt.Die Grundidee des Beweises sieht man im diskreten Fall:

i) =⇒ ii) folgt direkt aus der Definition der bedingten Wahrscheinlichkeit:

fϑ(x) = Ws({X = x} ∩ {V (X) = V (x)})= Ws({V (X) = V (x)})︸︷︷︸

=:gϑ(V (x))

·Ws({X = x}|{V (X) = V (x)})︸︷︷︸hangt nach i) nicht von ϑ ab

Fur die umgekehrte Richtung ist zu zeigen, dass

Wsϑ({X = x}|{V (X) = V (x)}) =Ws({X = x})

Ws({V (X) = V (x)})

nicht von ϑ abhangt. Wir verwenden dafur die Faktorisierung und kurzen den Faktorgϑ(V (x)):

Ws({V (X) = V (x)}) =∑

y|V (y)=V (x)

Ws(X = y)ii)= gϑ(V (x))

∑y|V (y)=V (x)

h(y)

Ws({X = x}) = gϑ(V (x)) · h(x)

Also hangt

Wsϑ({X = x}|{V (X) = V (x)}) =h(x)∑

y|V (y)=V (x) h(y)

nicht mehr von ϑ ab.

Folgerungeni) Der Maximum-Likelihood-Schatzer ϑML fur ϑ ist eine Funktion der suffizienten Sta-tistik: Da h(x) nicht von ϑ abhangt, genugt es zur Bestimmung von ϑML, gϑ(V (X)) zumaximieren.ii) Mit Hilfe des Faktorisierungssatzes kann man anhand der Dichte haufig leicht suffi-ziente Statistiken erschließen:

Beispiel: Gammaverteilung X1, . . . , Xn seien unabhangig und Gamma(α, λ)-verteilt,d.h. fur jedes X gemaß der Dichte

f(x) = c(α, λ) · xα−1e−λx, x ≥ 0.

Die gemeinsame Dichte von X1, . . . , Xn ist dann

fα,λ(x1, . . . , xn) = c′(α, λ) · (x1 · . . . · xn)α−1 · e−λ(x1+···+xn)

Man sieht sofort: (∏Xi,

∑Xi) ist suffizient fur (α, λ).

Mit Hilfe einer suffizienten Statistik kann man manchmal Schatzer mit kleineremmittleren quadratischen Fehler erzeugen, wenn man die bedingte Erwartung verwendet.Der Einfachheit halber betrachten wir hier nur den diskreten Fall.

58

Die bedingte Erwartung Y sei eine ZVe mit diskretem Wertebereich, X habe Dichteoder Gewichte f(x). Der bedingte Erwartungswert vonX, gegeben das Ereignis {Y = k},ist eine Zahl

E[X|Y = k] =1

Ws({Y = k})E[I{Y=k}X].

Fur jeden Wert k, den Y annehmen kann, ist E[X|Y = k] eine Zahl. Der Zufall kommtnun durch Y ins Spiel: Die bedingte Erwartung von X, gegeben Y , E[X|Y ], ist eineZufallsvariable der Gestalt

E[X|Y ] =∑k

I{Y=k}E[X|Y = k].

Beispiel Ein zweistufiges Zufallsexperiment:i) Y = 1/3 mit Wahrscheinlichkeit 1/2, und Y = 2/3 mit Wahrscheinlichkeit 1/2.ii) Gegeben Y = y, sei X die Anzahl der Erfolge beim n-maligen y-Munzwurf, alsoX ∼ Bi(n, y). Es gilt

E[X|Y = 1/3] = 1/3 · n, E[X|Y = 2/3] = 2/3 · n.

Damit ist die bedingte Erwartung E[X|Y ] vonX gegeben Y die Zufallsvariable E[X|Y ] =Y · n mit Werten in {n/3, 2n/3}, jeweils mit Wahrscheinlichkeit 1/2.

Der Satz von Rao-Blackwell Wieder beschranken wir uns auf den diskreten Fall,Analoges gilt fur den Dichtefall.Sei S = h(X) ein Schatzer fur ϑ. Sei V (X) eine suffiziente Statistik fur ϑ. Dann gilt furden Schatzer S∗ mit

S∗ := E[S|V (X)] (die ‘Rao-Blackwellisierung’ von S)

i) Die Erwartungswerte von S und S∗ sind gleich.ii) Der mittlere quadratische Fehler von S∗ ist hochstens so groß wie der von S.

Bew.i) Es gilt allgemein

E[E[X|Y ]] = E

[∑k

I{Y=k}E[X|Y = k]

]

=∑k

Ws({Y = k}) ·E[I{Y=k}X]

Ws({Y = k})

=∑k

E[I{Y=k}X] = E[X∑

I{Y=k}

]= EX,

also auchE[S∗] = E[E[S|V (X)]] = E[S].

ii) z.z.: E[(S∗ − ϑ)2] ≤ E[(S − ϑ)2]. Zunachst ist

S∗ − ϑ = E[S|V (X)]− ϑ = E[(S − ϑ)|V (X)]

E[(S∗ − ϑ)2] = E(E[(S − ϑ)|V (X)]2).

59

Um zu zeigen, dass E(E[(S − ϑ)|V (X)]2) ≤ E[(S − ϑ)2], beobachten wir:

E[(S − ϑ)|V (X)]2 =

(∑k

I{V (X)=k}E[S − ϑ|V (X) = k]

)2

=∑k

I{V (X)=k}E[S − ϑ|V (X) = k]2

≤∑k

I{V (X)=k}E[(S − ϑ)2|V (X) = k]

= E[(S − ϑ)2|V (X)]

und damitE[E[(S − ϑ)|V (X)]2] ≤ E[(S − ϑ)2].

�

Beispiel Seien X1, . . . , Xn unabhangige p-Munzwurfe. Wir betrachten die StatistikS := X1 ·X2 zur Schatzung von p2. Es gilt offensichtlich: EpS = p2, und wir wissen, dassdie Anzahl der Erfolge, K :=

∑Xi suffizient ist fur p. Also ist die Rao-Blackwellisierung

von S gegeben durch:S∗ := E[X1 ·X2|K]

MitE[X1 ·X2|K = k] = Ws(X1 = 1, X2 = 1|K = k) =

k

n· k − 1n− 1

folgt

S∗ =K(K − 1)n(n− 1)

.

Bemerkung Suffiziente Statistiken konnen offensichtlich nutzlich sein fur die effizienteBeschreibung eines Datensatzes sowie zur Erstellung von geeigneten Schatzern. Mansollte allerdings beachten, dass man damit die Daten stark reduziert und potentielleAbweichungen von den Modellannahmen nicht mehr erkennen kann.

Beispiel: X1, . . . , X100 unabhangige p-Munzwurfe. K :=∑Xi = 50. Wir schatzen

p = 1/2. Wenn aber alle Erfolge am Anfang der Serie auftraten und alle Misserfol-ge am Ende, haben wir durchaus Grund, an den Modellannahmen von unabhangigenMunzwurfen mit demselben Erfolgsparameter zu zweifeln. Dies ist nicht zu erkennen,wenn man die Daten auf die Anzahl K der Erfolge reduziert.

60

8 Likelihood-Quotienten-Tests

Die Likelihoodfunktion wird nicht nur eingesetzt, um Parameter zu schatzen, sondernauch um zu quantifizieren, wie gut Daten zu einer Hypothese (einer Annahme ubereinen bestimmten Parameter) passen, und um verschiedene Modelle / Hypothesen zuvergleichen.

Beispiel Wir betrachten zunachst ein klassisches Beispiel aus der Qualitatskontrolle:Eine Lieferung, bestehend aus n Einzelstucken, ist zu testen. Ein Anteil p der Ein-zelstucke ist schlecht, der Rest ist gut. Konsument und Produzent sind sich (z.B.) einig,dass p ≤ 0.03 annehmbar ist, aber p ≥ 0.15 nicht akzeptabel.Es werden 20 Stuck rein zufallig gezogen, die Anzahl X der schlechten Stucke wirdnotiert. Also X ∈ B = {0, 1, . . . , 20}. Man einigt sich auf:

A0 := {0, 1, 2} : Die Lieferung wird gekauft, falls X ∈ A0

A1 := B \A0 := {3, 4, . . . , 20} : Die Lieferung wird abgelehnt, falls X ∈ A1

Das Risiko des Produzenten und des Konsumenten haben wir bereits kennengelernt alsRisiko erster und zweiter Art:

Wsp≤0.03(X ∈ A1) ≤Wsp=0.03(X > 2) = 1− pbinom(2, 20, 0.03) ≈ 0.021.Wsp≥0.15(X ∈ A0) ≤Wsp=0.15(X ≤ 2) = pbinom(2, 20, 0.15) ≈ 0.4.

Der allgemeine Fall Gegeben sei eine ZVe X mit Werten in B, verteilt gemaß derDichte (oder der Gewichte) fϑ(x). Es sei ϑ ∈ Θ, wobei Θ die Menge der Parameterbezeichnet (im Beispiel Θ = [0, 1]). Seien Θ0,Θ1 zwei disjunkte Teilnehmen von Θ (imBsp. Θ0 = [0, 0.003], Θ1 = [0.15, 1]). Ein Test von Θ0 gegen Θ1 teilt den Wertebereichvon X in zwei disjunkte Teilmengen A0, A1 mit A0 ∪ A1 = B — die Annahmebereichefur die Hypothesen ϑ ∈ Θ0 bzw. ϑ ∈ Θ1. Man entscheide dann fur Θ0, falls X ∈ A0,sonst fur Θ1. Wie gehabt:

Risiko 1. Art: supϑ∈Θ0

Wsϑ(X /∈ A0) =: α

Risiko 2. Art: supϑ∈Θ1

Wsϑ(X /∈ A1) =: β.

Die Abbildungϑ→Wsϑ(X ∈ A0)

heißt Operationscharakteristik des Tests. Fur ϑ ∈ Θ0 soll sie moglichst groß, sonstmoglichst klein sein.

Gesucht ist nun ein Test, der fur vorgegebenes Risiko 1. Art (d.h. fur vorgegebenesSignifikanzniveau) das kleinste Risiko 2. Art besitzt (also die großte Testmacht). ImSpezialfall ist das einfach:

61

8.1 Einfache Hypothese und einfache Alternative

Seien Θ0 = {ϑ0}, Θ1 = {ϑ1}. Seien f0(x), f1(x) die Verteilungsdichten (bzw. -gewichte)von X unter ϑ0 bzw. ϑ1. Wir betrachten den Likelihood-Quotienten

q(x) :=f0(x)f1(x)

.

Sind die Beobachtungen unter ϑ0 wahrscheinlicher als unter ϑ1, so ist q(x) > 1, sonstq(x) ≤ 1. Je großer das Verhaltnis dieser Gewichte wird, desto mehr sprechen die Datenfur ϑ0 (und umgekehrt).

Das Lemma von Neyman und Pearson Sei c > 0 fest und so gewahlt, dass

A(c) :={x ∈ B

∣∣∣∣f0(x)f1(x)

≥ c}

ein Annahmebereich fur ϑ0 zum Niveau 1 − α ist und sein Niveau voll ausschopft, sodass

Wsϑ0(X /∈ A(c)) = α.

Dann ist der Test mit Annahmebereich A0 = A(c) unter allen Tests von {ϑ0} gegen{ϑ1} mit Signifikanz α der Machtigste.

Bew.:Sei A ein weiterer Annahmebereich fur {ϑ0} zum Niveau 1− α. Dann ist zu zeigen:

Wsϑ1(X ∈ A) ≥Wsϑ1(X ∈ A(c)).

Wir betrachten zunachst

Wsϑ1(X ∈ A) =∫Af1(x)dx =

∫A(c)

f1(x)dx+∫A\A(c)

f1(x)dx−∫A(c)\A

f1(x)dx.

Auf der Menge A \A(c) ist f1 > f0/c, und auf A(c) ist f1 ≤ f0/c. Daher:

Wsϑ1(X ∈ A) ≥∫A(c)

f1(x)dx+1c

∫A\A(c)

f0(x)dx− 1c

∫A(c)\A

f0(x)dx

=∫A(c)

f1(x)dx+1c

∫Af0(x)dx︸︷︷︸

=Wsϑ0(X∈A)≥1−α

−∫A(c)

f0(x)dx︸︷︷︸=Wsϑ0

(X∈A(c))=1−α

≥∫A(c)

f1(x)dx = Wsϑ1(X ∈ A(c)).

�

62

8.2 Generalisierte Likelihood-Quotiententests

Tests nach Neyman und Pearson testen eine einfache Hypothese (ein moglicher Wertfur ϑ) gegen eine einfache Alternative. Daher war es relativ leicht, ihre Optimalitat zuzeigen. Likelihood-Quotiententests sind allerdings viel breiter anwendbar, und obwohlnicht zwangslaufig optimal, so funktionieren sie haufig gut und spielen daher beim Testeneine ahnliche Rolle wie das Maximum-Likelihood-Prinzip in der Schatztheorie.

Wir betrachten unser Beispiel vom Anfang: Seien Θ0,Θ1 zwei disjunkte Teilmengenvon Θ (d.h. nicht mehr zwangslaufig einelementig, wie bei Neyman-Pearson), mit Θ0 ∪Θ1 = Θ. Dann ist der verallgemeinerte Likelihood-Quotient definiert als

Q∗(x) :=supϑ∈Θ0

fϑ(x)supϑ∈Θ1

fϑ(x).

Aus technischen Grunden bevorzugt man meistens die folgende Große:

Q(x) :=supϑ∈Θ0

fϑ(x)supϑ∈Θ fϑ(x)

,

wobei Q(x) = min(Q∗(x), 1). (Denn ist der ML-Schatzer fur ϑ in Θ0, so ist Q(x) = 1und Q∗(x) ≥ 1, ist der ML-Schatzer in Θ1, so ist Q(x) = Q∗(x).)

Große Werte von Q sprechen fur H0 : ϑ ∈ Θ0. Daher ist analog zu Neyman-Pearsonder Annahmebereich A0 fur H0 festgelegt durch alle Beobachtungen, fur die Q großerist als eine Schranke c: A0 := {x|Q(x) ≥ c}, fur ein c ∈ R.

Beispiel: Normalverteilung Seien X1, . . . , Xn unabhangig und identisch N (µ, σ2)-verteilt, mit bekanntem σ2. Sei H0 : µ = µ0 und H1 : µ 6= µ0, d.h. Θ0 = {µ0},Θ1 = R \ {µ0}. Der Zahler von Q(X) ist damit einfach gegeben durch

1(√

2πσ)ne−

12σ2

P(Xi−µ0)2 .

Der Nenner ist der Wert von f am ML-Schatzer X von µ :

1(√

2πσ)ne−

12σ2

P(Xi−X)2 .

Damit ist der Likelihood-Quotient gegeben durch

Q(X) = exp(− 1

2σ2

(∑(Xi − µ0)2 −

∑(Xi − X)2

))Q wird klein, wenn −2 logQ groß wird:

−2 logQ(X) =1σ2

(∑(Xi − µ0)2 −

∑(Xi − X)2

)=n(X − µ0)2

σ2.

Hier treffen wir einen alten Bekannten, denn −2 logQ ist gerade das Quadrat der z-Statistik zum Test der Hypothese H0 : µ = µ0. Z ∼ N (0, 1), also −2 logQ(X) ∼ χ2(1).

63

Asymptotische Verteilung von Q Wenn die Verteilung von Q unter H0 bekannt ist,ist es leicht, einen Annahme- und Ablehnungsbereich fur H0 zu definieren und damitH0 gegen ihre Alternative zu testen. Die Verteilung von Q ist nicht immer bekannt,aber unter bestimmten Regularitatsvoraussetzungen an die Dichte von X ist −2 logQasymptotisch χ2-verteilt mit (dim(Θ) − dim(Θ0)) Freiheitsgraden. Im obigen Beispielsogar exakt: −2 logQ ∼ χ2(1) = χ2(1− 0) = χ2(dim(Θ)− dim(Θ0)).

Beispiel: Multinomialverteilung Auch hier ist der Likelihood-Quotiententest (zu-mindest asymptotisch) aquivalent zu einem haufig eingesetzten Test, den wir im nachstenAbschnitt diskutieren werden: dem χ2-Test zum Vergleich beobachteter und erwarteterHaufigkeiten.

Y1, . . . , Yn nehmen Werte in k Kategorien an. X1, . . . , Xk seien die beobachtetenHaufigkeiten der k Kategorien. Dann ist X = (X1, . . . , Xk) multinomialverteilt mit Pa-rametern n und (p1, . . . , pk). Die Menge der moglichen Parameter ist

Θ = {(p1, . . . , pk)|pi ≥ 0,∑

pi = 1}.

Wir untersuchen wieder eine einfache Hypothese:

Θ0 := {(π1, . . . , πk)}, mit festen Gewichten π1, . . . , πk,∑

πi = 1.

Fur den Nenner vonQ brauchen wir die ML-Schatzer der Zellbesetzungen (p1, . . . , pk).Dies sind gerade die relativen Haufigkeiten, pi = Xi/n (Ubung).

Der Likelihood-Quotient ist damit

Q(X) =n!

X1!···Xk!πX11 · · ·π

Xkk

n!X1!···Xk! p

X11 · · · p

Xkk

=∏(

πipi

)Xi.

Also:

−2 logQ(X) = −2∑

Xi log(πipi

)= 2n

∑pi log

(piπi

).

Auch diese Große ist entsprechend obiger Bemerkung approximativ χ2(dim(Θ)−dim(Θ0)) =χ2(k−1)−verteilt (dim Θ = k−1 wegen

∑pi = 1). Außerdem ist −2 logQ asymptotisch

gleich der Pearsonschen Teststatistik X fur die Hypothese Hπ : p = (π1, . . . , πk):

X :=∑ (Xi − nπi)2

nπi=:∑ (Bi − Ei)2

Ei,

wobei Bi := Xi die beobachteten Haufigkeiten, und Ei := nπi die erwarteten Haufigkei-ten in den k Kategorien bezeichnen.

Die asymptotische Gleichheit von X und−2 logQ sieht man an einer Taylor-Entwicklung:Fur x nahe x0 ist

g(x) = x log(x

x0

)approximierbar durch

g(x) ≈ g(x0) + (x− x0)g′(x0) +(x− x0)2

2g′′(x0) = (x− x0) +

(x− x0)2

2x0.

64

Damit ist (fur p = (p1, . . . , pk) nahe π = (π1, . . . , πk), d.h. fur genugend großes n unterder Nullhypothese p = π)

−2 logQ = 2n∑

pi log(piπi

)≈ 2n

∑((pi − πi) +

(πi − pi)2

2πi

)= 2n

∑(pi − πi) + n

∑ (πi − pi)2

πi=∑ (npi − nπi)2

nπi= X .

9 Der χ2-Test

Wir betrachten die Pearsonsche Teststatistik fur feste hypothetische Gewichte:

Satz (Pearson, 1900)X(n) = (X(n)

1 , . . . , X(n)k ) sei multinomial (n;π1, . . . , πk)-verteilt. Dann gilt:

X :=k∑i=1

(X(n)i − nπi)2

nπi−→ χ2(k − 1) fur n→∞ in Verteilung.

Die erste Plausibilitatsuberlegung hierzu ist die folgende: Fur jedes feste i ist X(n)i ∼

Bi(n, πi), also gilt nach dem Satz von de Moivre-Laplace:(X

(n)i − nπi√nπi(1− πi)

)2

−→ χ2(1) fur n→∞ in Verteilung.

Waren alle X(n)i unabhangig, so wurde gelten

∑i

(X

(n)i −nπi√nπi(1−πi)

)2

→ χ2(k). Ein Freiheits-

grad geht verloren durch die Abhangigkeit derX(n)i (denn ihre Summe ist deterministisch

n). Die passende Varianz im Nenner sieht man beim Ubergang zur Poissonverteilung:

Fishers heuristischer Beweis Wir poissonisieren: Seien N (n)1 , . . . , N

(n)k unabhangig

mit N (n)i ∼ Po(nπi). Dann ist die Summe der Ni nicht mehr deterministisch n, sondern

poissonverteilt mit Parameter n:

N (n) :=∑

N(n)i ∼ Po(n).

Die ZVe N (n)i ∼ Po(nπi) hat Varianz nπi, und daher folgt mit dem zentralen Grenz-

wertsatz:N

(n)i − nπi√nπi

−→ N (0, 1) fur n→∞ in Verteilung.

Damit konvergiert der Vektor Gn mit Eintragen (N(n)i −nπ√nπi

)i in Verteilung gegen eine

standardnormalverteilte Zufallsvariable im Rk:

Gn :=

(N

(n)i − nπi√nπi

)i=1,...,k

−→ Z fur n→∞ in Verteilung.

65

Aus dem Satz von Fisher wissen wir, dass die Projektionen von Z auf orthogonale Un-terraume U1,U2 von Rk unabhangig sind und ihre Langenquadrate χ2-verteilt mit Frei-heitsgraden dim(U1) und dim(U2). Wir zerlegen Gn, was ja in Verteilung gegen Z konver-giert, durch orthogonale Projektion (und projizieren damit sozusagen den zusatzlichenFreiheitsgrad der zufalligen Gesamtzahl N (n) weg). Projektion auf den von

e :=

√π1...√πk

aufgespannten Unterraum liefert

eTGn = 〈e,Gn〉 =1√n

(N (n) − n) −→ eTZ fur n→∞ in Verteilung,

wobei nach dem Satz von Fisher eTZ ∼ N (0, 1).Wir zerlegen Gn damit in die Projektion PeGn auf e und die Projektion Pe⊥Gn auf

das orthogonale Komplement von e:

Gn = PeGn + Pe⊥Gn,

PeGn = (eTGn)e =(

1√n

(N (n)√πi − n√πi))i=1,...,k

Pe⊥Gn = Gn − PeGn =1√n

(N

(n)i√πi−N (n)√πi

)i=1,...,k

Es gilt Pe⊥Gn → Pe⊥Z in Verteilung, und damit

||Pe⊥Gn||2 =∑ (N (n)

i −N (n)πi)2

nπi−→ χ2(k − 1).

Außerdem ist ||Pe⊥Gn||2(→ ||Pe⊥Z||2) asymptotisch unabhangig von ||PeGn||2(→ ||PeZ||2) =1n(N (n) − n)2 und damit von N (n).

Abschließend beobachten wir die Beziehung zwischen der Poisson- und der Multino-mialverteilung: Bedingt unter {N (n) = n} ist (N (n)

1 , . . . , N(n)k ) multinomial(n, π1, . . . , πk)-

verteilt:

Ws({N (n)1 = y1, . . . , N

(n)k = yk}|{N (n) = n})

=Ws({N (n)

1 = y1, . . . , N(n)k = yk} ∩ {N (n) = n})

Ws({N (n) = n})

=e−nπ1 (nπ1)y1

y1! · · · e−nπk (nπk)yk

yk!

e−n nn

n!

=(

ny1 · · · yk

)πy11 · · ·π

ykk .

Fur einen rigorosen Beweis brauchen wir Wissen uber Normalverteilungen auf Rk:

66

9.1 Normalverteilungen auf Rk

Kovarianzmatrizen Y = (Y1, . . . , Yk)T sei Rk-wertige Zufallsvariable mit Erwar-tungswert m := EY . Sei Y 0 := Y −m. Dann ist die Matrix C mit Eintragen

Cij := Cov(Yi, Yj) = E[(Yi −mi)(Yj −mj)] = E[Y 0i Y

0j ]

die Kovarianzmatrix von Y .

Jede Kovarianzmatrix ist symmetrisch und positiv semidefinit.Denn fur a = (a1, . . . , ak)T ∈ Rk ist

aTCa =∑i,j

aiCov(Yi, Yj)aj = Cov(∑

aiYi,∑

ajYj)

= Var(∑

aiYi)≥ 0.

Jede symm., positiv semidefinite Matrix tritt als Kovarianzmatrix auf.C sei k× k-Matrix, symmetrisch und positiv semidefinit. Dann gibt es Eigenwerte σ2

1 ≥σ2

2 ≥ . . . ≥ σ2k ≥ 0 und eine zugehorige ONB aus Eigenvektoren (Beweis: Lineare

Algebra.) Das heißt, es gibt eine Darstellung von C der Gestalt:

C =k∑i=1

σ2i eie

Ti .

Denn in der Tat ist Cei = σ2i ei = (

∑kj=1 σ

2j eje

Tj )ei.

DefinitionY heißt normalverteilt auf Rk mit Mittelwert m und Kovarianzmatrix C (Y ∼ N (m,C))⇐⇒ Y hat Erwartungswert m und Kovarianzmatrix C, und ∀a ∈ Rk : aTY ist normal-verteilt auf R.

SatzSei C =

∑ki=1 σ

2i eie

Ti . Dann gilt:

N (0, C) = L( k∑i=1

σiZiei),

wobei Z1, . . . , Zk unabhangig und N (0, 1)-verteilt auf R. Denn mit Y :=∑σiZiei ist

E[Y Y T ] = E[(∑i

σiZiei)(∑j

σjZjej)] =∑i

σ2i eie

Ti .

Beispiel C =∑r

i=1 eieTi mit r < k, ei orthonormal. Dann ist C eine Projektionsmatrix

auf einen r-dimensionalen Teilraum von Rk. Dann ist N (0, C) die Standardnormalver-teilung auf dem von e1, . . . , er aufgespannten Teilraum von Rk.

67

Wie im Eindimensionalen gilt auch im Rk der (nun multivariate) zentrale Grenz-wertsatz: Seien Y1, . . . , Yn u.i.v. im Rk mit Mittelwert 0 und Kovarianzmatrix C. Danngilt:

1√n

n∑i=1

Yi −→ N (0, C) fur n→∞ in Verteilung.

(Beweis Vorlesung Hohere Stochastik, Fr. 14.12.07, 10-12 Uhr).Beweis (Satz von Pearson):Seien Y1, . . . , Yn u.i.v. mit Werten in {1, . . . , k} und Verteilungsgewichten (πi)i=1,...,k.Seien fur i = 1, . . . , n

Vi :=

1√π1I{Yi=1}...

1√πkI{Yi=k}

.

Dann sind alle Vi unabhangig identisch verteilt mit Erwartungswert e = (√π1, . . . ,

√πk)T

und Kovarianzmatrix

Cij = E[(

1√πiI{Y=i} −

√πi

)·(

1√πjI{Y=j} −

√πj

)]= δij −

√πiπj ,

also C = I − eeT . Damit folgt

1√n

(∑Vi − ne

)=

1√n

1√π1X1 − n

√π1

...1√πkXk − n

√πk

→ N (0, I−eeT ) fur n→∞ in Verteilung.

Also gilt ∑(Xi − nπi√

nπi

)2

=∣∣∣∣ 1√

n

(∑Vi − ne

) ∣∣∣∣2 → ||Z||2mit Z ∼ N (0, I − eeT ). Laut obigem Beispiel ist Z damit standardnormalverteilt aufRk−1, denn I − eeT ist Projektionsmatrix auf einen (k − 1)-dimensionalen Unterraum.Damit ist die Pearsonsche Statistik X approximativ χ2(k − 1)-verteilt. �

9.2 Der χ2-Test

Satz (Pearson, Fisher)(X1, . . . , Xk) seien beobachtete Haufigkeiten und multinomial (n, π)-verteilt. Gegebensei eine Hypothese: π ∈ F , wobei F eine r-dimensionale Hyperflache des Simplex

Θ := {(p1, . . . , pk)|∑

pi = 1, pi ≥ 0}

ist. pF sei der ML-Schatzer fur π unter der Hypothese π ∈ F . Dann ist die Teststatistik

T :=k∑i=1

(Xi − npFi )2

npFi

unter π ∈ F fur große n annahernd χ2(k − 1− r)-verteilt.

68

Beispiel: k = 4. Y1, . . . , Yn fallen in 4 mogliche Kategorien, mit

Θ ={(

p11 p12

p21 p22

) ∣∣∣∣∑ pij = 1, pij ≥ 0}

Offensichtlich gilt dim(Θ) = 3. Sei nun die Hypothese

F :={(

αβ α(1− β)(1− α)β (1− α)(1− β)

) ∣∣∣∣0 ≤ α ≤ 1, 0 ≤ β ≤ 1},

mit dim(F ) = 2. Was sind die ML-Schatzer fur α und β fur gegebene Haufigkeitenx11, x12, x21 und x22? Wir maximieren die Multinomialgewichte:

fαβ(x) =(αβ)x11 · (α(1− β))x12 · ((1− α)β)x21 · ((1− α)(1− β))x22 · c=αx1+(1− α)x2+βx+1(1− β)x+2 · c,

mit xi+ :=∑

j xij und x+i :=∑

j xji. Durch Ableiten der log-likelihood sieht man:

α =x1+

nβ =

x+1

n.

Der χ2−Test auf UnabhangigkeitDas obige Beispiel mit k = 4 ist ein Spezialfall fur den χ2−Test auf Unabhangigkeit.Gegeben sind dabei zwei Merkmale A und B in I bzw. J Abstufungen, also insgesamtk = I ·J Kategorien. Wir beobachten die Haufigkeiten xij in den Kategorien (ij). pij seidie Wahrscheinlichkeit, dass ein zufallig gezogenes Individuum in Kategorie (ij) liegt.Y1, . . . , Yn seien u.i.v. mit Gewichten pij . Die Hypothese ist:

pij = αiβj ,

wobei αi die Wahrscheinlichkeit bezeichnet, in Merkmal A in Kategorie i zu liegen,und βj entsprechend die Wahrscheinlichkeit fur Kategorie j in Merkmal B, mit

∑αi =∑

βj = 1. Die Hypothese sagt also, dass die beiden Merkmale in der Population un-abhangig sind, mit dim(F ) = (I − 1) + (J − 1) = I + J − 2. Die ML-Schatzer ermitteltman analog zu oben:

αi :=xi+n

βj =x+j

n,

mit xi+ und x+j wie oben. Dann folgt mit dem Satz von Pearson und Fisher: Ist(X11, X12, . . . , XIJ) multinomial (n, α1β1, α1β2, . . . , αIβJ)-verteilt, so ist

∑i,j

(Xij − nαiβj)2

nαiβj=∑(

Xij − Xi+X+j

n

)2Xi+X+j

n

fur große n approximativ χ2(IJ − 1− (I + J − 2)) = χ2((I − 1)(J − 1))-verteilt.

69

Statistik - uni-frankfurt.de · 1 Beschreibende Statistik Ein wesentlicher Teil der Statistik besch...

Documents

Transcript of Statistik - uni-frankfurt.de · 1 Beschreibende Statistik Ein wesentlicher Teil der Statistik besch...