Generalisierte Lineare Mo delle mit zuf
alligen
Eekten und variierenden Ko eÆzienten
vorgelegt von
Diplom-Wirtschaftsmathematikerin
Claudia Beate Funck-H
usges
Von der Fakult
at IV - Elektrotechnik und Informatik
der Technischen Universit
at Berlin
zur Erlangung des akademischen Grades
Doktorin der Naturwissenschaften
- Dr.rer.nat. -
genehmigte Dissertation
Promotionsausschuss:
Vorsitzender: Prof. Dr.-Ing. Adam Wolisz
Berichter: Prof. Dr. Gerhard Tutz
Berichter: Prof. Dr. UlrichKockelkorn
Tag der wissenschaftlichen Aussprache: 03.12.2001
Berlin 2001
D83
Inhaltsverzeichnis
Einleitung 1
1 Das Generalisierte Lineare Mo dell 3
1.1 Grundlagen des Generalisierten Linearen Mo dells . . . . . . . 4
1.1.1 Mo delldenitionen . . . . . . . . . . . . . . . . . . . . 4
1.1.2 Sch
atzung der Parameter . . . . . . . . . . . . . . . . . 11
2 Mo dellerweiterungen im GLM 13
2.1 Mo dellans
atze mit zuf
alligen Eekten . . . . . . . . . . . . . . 16
2.1.1 EM-Algorithmus mit Gau-Hermite Integrationstechnik 19
2.1.2 EM-Typ e Algorithmus . . . . . . . . . . . . . . . . . . 25
2.1.3 Mo dellansatz der nichtparametrischen Maximum
Likeliho o d . . . . . . . . . . . . . . . . . . . . . . . . . 28
2.2 Mo dellans
atze mit (zeit-)variierenden Ko eÆzienten . . . . . . 33
2.2.1 Parametrische und semi-parametrische Ans
atze . . . . 34
2.2.2 Lokaler Likeliho o dansatz . . . . . . . . . . . . . . . . . 36
3Zuf
allige Eekte und variierende Ko eÆzienten 40
3.1 Mo dellans
atze und Parametersch
atzungen . . . . . . . . . . . 42
3.1.1 Volles Mo dell, parametrische und semi-
parametrische Ans
atze . . . . . . . . . . . . . . . . . . 42
3.1.2 Lokaler Likeliho o dansatz mit zuf
alligen Eekten . . . . 43
i
3.2 Varianz der Parametersch
atzungen . . . . . . . . . . . . . . . 48
3.2.1 Varianzsch
atzungen b ei Mo dellen mit zuf
alligen Eekten 48
4 Simulationsb eispiele 54
4.1 Sch
atzung der festen Parameter . . . . . . . . . . . . . . . . . 56
4.1.1 Das volle Mo dell mit zuf
alligen Eekten . . . . . . . . 56
4.1.2 Lokaler Likeliho o dansatz mit zuf
alligen Eekten und
Gau-Hermite Integrationstechnik . . . . . . . . . . . . 60
4.1.3 Lokaler Likeliho o dansatz mit zuf
alligen Eekten und
EM-Typ e Algorithmus . . . . . . . . . . . . . . . . . . 64
4.2 Bayes-Sch
atzung der zuf
alligen Eekte . . . . . . . . . . . . . 67
5 Anwendungsb eispiele 70
5.1 Poissonmo delle . . . . . . . . . . . . . . . . . . . . . . . . . . 71
5.1.1 Emesis Datensatz . . . . . . . . . . . . . . . . . . . . . 71
5.1.2 Epileptic Datensatz . . . . . . . . . . . . . . . . . . . . 77
5.2 Verweildauermo delle . . . . . . . . . . . . . . . . . . . . . . . 80
5.2.1 Dauer der Arb eitslosigkeit . . . . . . . . . . . . . . . . 80
5.2.2 Magenkrebsstudie . . . . . . . . . . . . . . . . . . . . . 82
6 Anhang 87
6.1 EM-Algorithmus mit direkter Gau-Hermite Integrationstechnik 87
6.1.1 Numerische Integration mit Gau-Hermite . . . . . . . 87
6.1.2 Sch
atzung der festen Parameter . . . . . . . . . . . . . 89
6.2 EM-Typ e Algorithmus . . . . . . . . . . . . . . . . . . . . . . 97
6.2.1 Blo ckdiagonalgestalt der Fishermatrix . . . . . . . . . 97
Zusammenfassung 99
Literaturverzeichnis 102
ii
Einleitung
Generalisierte Lineare Mo delle (GLM) sind ein vielverwendetes Instrumenta-
rium f
ur die Durchf
uhrung von Zusammenhangsanalysen, die b eispielsweise
bei den Auswertungen von medizinischen o der sozial
okonomischen Studien
verwendet werden.
Ihre groe Verbreitung liegt zum einen daran, da in dieser Mo dellklasse
sowohl das Lineare Mo dell als auch Sp ezialf
alle f
ur z. B. diskrete, kategoriale
oder p oissonverteilte Daten enthalten sind und somit sehr unterschiedliche
Datenmo delle analysiert werden k
onnen.
Andererseits sind Generalisierte Lineare Mo delle in den verschiedensten Sta-
tistik-Software Programmen implementiert und sind dadurch sehr einfach
vom Anwender abrufbar.
Das grundlegende Generalisierte Lineare Mo dell gehtvon unabh
angigen Be-
obachtungseinheiten aus. Jedo ch liegen in vielen groangelegten Studien
Grupp enstrukturen vor, so da sp ezielle Mo dellerweiterungen angewendet
werden m
ussen. Diese zu ber
ucksichtigenden Grupp enstrukturen k
onnen
durch Mewiederholungen an den gleichen Beobachtungseinheiten o der auch
durch eine famili
are, regionale o der administrative Zugeh
origkeit der Beob-
achtungseinheiten gegeb en sein.
Zus
atzlich zu der Grupp enstruktur mu h
aug auchnocheinm
oglicher zeit-
licher Einu ber
ucksichtigt werden, da die Individuen oder Beobachtungs-
1
einheiten in diesen Studien meist
ub er einen l
angeren Zeitraum b eobachtet
werden.
Das Ziel dieser Arb eit ist es, die grupp en- und/o der individuensp ezische Da-
tenstruktur durch die Aufnahme von zuf
alligen Eekten und gleichzeitig den
zeitlichen Einusses durchBetrachtung von zeitvariierenden Ko eÆzienten zu
ber
ucksichtigen.
Nach der Einleitung in die Thematik dieser Arb eit wird im ersten Kapitel das
grundlegende Generalisierte Lineare Mo dell vorgestellt. Auf die Mo dellerwei-
terungen mit zeitvariierenden Ko eÆzienten und mit zuf
alligen Eekten wird
im zweiten Kapitel eingegangen. Diese b eiden Mo dellans
atze werden im drit-
ten Kapitel in gemeinsame Mo delle aufgenommen. Im vierten Kapitel wird
in den Simulationsb eispielen neb en der G
ute f
ur die Parametersch
atzungen
selbst auch die G
ute f
ur die Sch
atzung der Varianzen der Parameter be-
trachtet. Aus der Vielzahl der Anwendungsm
oglichkeiten werden im f
unften
Kapitel vier Anwendungsb eispiele herausgestellt. Das sechste Kapitel enth
alt
einige technische Details f
ur die im zweiten Kapitel b eschrieb enen Verfahren.
In der nachfolgenden Zusammenfassung erfolgt eine, diese Arb eit abschlie-
ende, Beurteilung der vorgestellten Verfahren.
2
Kapitel 1
Das Generalisierte Lineare
Mo dell
Das Generalisierte Lineare Mo dell (GLM), welches von Nelder & Wedderburn
(1972) vorgeschlagen wurde, ist eine Erweiterung des klassischen Regres-
sionsansatzes im Linearen Mo dell.
Die im folgenden vorgestellten Mo dellannahmen im Linearen Mo dell sind f
ur
viele metrische Zielgr
oen gut motiviert, falls ein linearer Zusammenhang
zwischen der Resp onsevariablen und den erkl
arenden Variablen, die auch als
Regressoren bzw. Kovariablen b ezeichnet werden, b esteht.
Liegt die Zielvariable diskret, b eispielsweise in Form von binomialen o der
multinomialen Zielgr
oen o der als Z
ahldaten vor, so ist das klassische lineare
Mo dell meist nicht mehr passend und es sollte b esser die Erweiterung durch
das Generalisierte Lineare Mo delle angewendet werden.
3
1.1 Grundlagen des Generalisierten Linearen
Mo dells
Da das Generalisierte Lineare Mo dell eine Erweiterung des Linearen Mo dells
ist, werden wir zun
achst kurz auf diesen Regressionsansatz eingehen und im
Anschlu die Mo delldenition des Generalisierten Linearen Mo dells vorstel-
len.
1.1.1 Mo delldenitionen
Das Lineare Mo dell
F
ur den univariaten Fall der linearen Regression wird vorausgesetzt, da
Beobachtungen (
y
i
;x
i
) mit
i
= 1
;:::;N
vorliegen. Bei dem Mo dellansatz
wird zwischen der eindimensionalen Resp onsevariablen
y
i
und dem Vektor
x
0
i
= (
x
i
1
;:::;x
im
)der erkl
arenden Variablen ein linearer Zusammenhang
der Form
y
i
=
z
0
i
+
i
f
ur
i
=1
;:::;N
angenommen.
Dab ei ist
0
= (
1
;:::;
p
) der Vektor des
p
-dimensionalen unbekannten
Parameter und
z
i
ein
p
-dimensionaler Designvektor, der durch eine geeig-
nete Funktion aus dem Kovariablenvektor
x
i
gebildet wird; b eispielsweise
durch
z
0
i
=(1
;x
0
i
). Weiter ist
i
eine Fehler- bzw. St
orvariable mit dem Er-
wartungswert
E
(
i
) = 0. Es wird vorausgesetzt, da
i
normalverteilt mit
i
N
(0
;
2
) ist und
y
i
b ei gegeb en
x
i
b edingt unabh
angig ist. Damit folgt
b ei gegeb enem
x
i
, da
y
i
N
(
i
;
2
)mit
i
=
E
(
y
i
j
x
i
) gilt.
Im linearen Mo dell ist
i
=
z
0
i
f
ur
i
=1
;:::;N
,das b edeutet, der bedingte
Erwartungswert ist gleich dem linearen Pr
adiktor
z
0
i
.
4
Das Generalisierte Lineare Mo dell
Im Generalisierten Linearen Mo dell wird im Gegensatz zum Linearen Mo-
dell die Verteilungsannahme f
ur die Zielgr
oe verallgemeinert indem voraus-
gesetzt wird, da die Dichte von
y
i
aus einer Exp onentialfamilie ist. Bei-
spiele f
ur die Exp onentialfamilie sind die Normal-, Binomial-, Multinomial-,
Poisson- und Gammaverteilung. Die Zielgr
oe
y
i
kann mehrdimensional sein,
d. h.
y
0
i
= (
y
i
1
;:::;y
iq
). Der Zusammenhang zwischen dem Erwartungs-
wert
E
(
y
i
j
x
i
)=
i
=(
i
1
;:::;
iq
)
0
und dem linearen Pr
adiktor
Z
i
mit der
Designmatrix
Z
i
und dem unbekannten Parametervektor
0
= (
1
;:::;
p
)
wird durch eine mehrdimensionale Resp onsefunktion
h
mit
i
=
h
(
Z
i
) her-
gestellt.
Das Generalisierte Lineare Mo dell, das b eispielsweise in McCullagh & Nelder
(1989) und Fahrmeir & Tutz (1994) vorgestellt ist, wird durchdie Vertei-
lungsannahme, die Strukturannahme und dem Aufbau der Designmatrix der
Kovariablen deniert.
Bei der Verteilungsannahme wird vorausgesetzt, da die
y
i
, gegeb en
x
i
,be-
dingt unabh
angig sind und ihre Verteilung aus einer Exp onentialfamilie
f
(
y
i
j
i
;;w
i
)=exp
(
y
0
i
i
b
(
i
))
w
i
c
(
y
i
;;w
i
)
:
stammt.
Durch die Strukturannahme wird der Erwartungswert
i
=(
i
1
;:::;
iq
)
0
=
E
(
y
i
j
x
i
) mit dem linearen Pr
adiktor
i
= (
i
1
;:::;
iq
)
0
=
Z
i
durch eine
invertierbare, zweimal stetig dierenzierbare Funktion
h
:
R
q
!
R
q
mit
i
=
h
(
i
)=
h
(
Z
i
)
verbunden.
5
Die Designmatrix
Z
i
2
R
q
p
wird durcheine geeignete Abbildung aus den
Kovariablen
x
i
=(
x
i
1
;:::;x
im
)
0
gebildet. Ist die Auspr
agung der Kovariablen
kategorial, so wird eine Ko dierung der Variablen vorgenommen.
Bemerkungen zur Verteilungsannahme
Die Funktionen
b
i
(
) und
c
i
(
) b estimmen den genauen Typ der Exp o-
nentialfamilie, d. h. ob b eispielsweise die Normal-, Binomial-, Multi-
nomial-, Poisson- o der Gammaverteilung vorliegt.
Der nat
urliche Parameter
i
ist eine Funktion des Erwartungswertes
i
,
d. h.
i
=
(
i
).
Der Erwartungswert ist durch die Wahl der Exp onentialfamilie mit
=
@b
(
)
@
und die Kovarianz mit
cov
(
y
i
j
x
i
)=
@
2
b
(
)
@ @
0
w
i
b estimmt.
Der Disp ersionsparameter
ist von
i
unabh
angig und in der Regel
unbekannt. Im Linearen Mo dell gilt
=
2
. Soweit nichts anderes
angegeb en ist, b etrachten wir im folgenden
als fest und bekannt mit
1.
Liegen die Daten in gruppierter Form vor und ist
N
i
die Anzahl der Be-
obachtungen in Grupp e
i
, so sind die Gewichte
w
i
=
N
i
f
ur
i
=1
;:::;N
falls die Zielvariable
y
i
den Grupp endurchschnitt b eschreibt, und es ist
w
i
=
1
N
i
, falls
y
i
die Summe der Auspr
agungen in den Beobachtungs-
grupp en ist. Sind die Daten nicht gruppiert, so gilt
w
i
1.
Bemerkungen zur Strukturannahme
Der Vektor
= (
1
;:::;
p
)
0
ist der unbekannte zu sch
atzende Para-
metervektor.
6
Die Funktion
h
(
i
) = (
h
1
(
i
1
)
;:::;h
q
(
iq
))
0
heit Resp onsefunktion.
Die zu
h
inverse Funktion
g
=
h
1
heit Linkfunktion mit
g
(
i
) =
i
=
Z
i
.
Der Werteb ereich der Resp onsefunktion
h
wurde in der Denition all-
gemein mit
R
q
deniert. Durch die Wahl der Linkfunktion wird der
Werteb ereich ab er meist geeignet eingeschr
ankt, z. B. f
ur die Multino-
mialverteilung auf (0
;
1)
(0
;
1).
W
ahlt man bei der Verteilungsannahme die Normalverteilung und f
ur die
Linkfunktion die identische Abbildung, so erh
alt man das Lineare Mo dell als
Sp ezialfall eines Generalisierten Linearen Mo dells.
Sp ezielle Datenmo delle
F
ur jede Exp onentialfamilie gibt es eine nat
urliche Linkfunktion, mit
g
(
)
(
), d. h. es gilt
Z
oder in Worten ausgedr
uckt, der nat
urliche Pa-
rameter ist gleichdem linearen Pr
adiktor. In dieser Arb eit werden wir Mo-
dellans
atze mit der nat
urliche Linkfunktion verwenden, die wir im folgenden
kurz vorstellen werden.
Poissonmo dell
Im Log-linearen Poissonmo dell wird
y
P
(
) angenom-
men und f
ur den unbekannten Parameter
wird die Mo dellannahme
=
E
(
y
j
x
) = exp(
z
0
) getroen, d. h. die Resp onsefunktion in diesem Mo dell
ist die Exp onentialfunktion. Der Designvektor
z
enth
alt den Regressoren-
vektor
x
und gegeb enenfalls auch die1f
ur den Interzept. Im Poissonmo dell
gilt f
ur die Varianz
var
(
y
j
x
)=
.
7
Binomiales Logit-Mo dell
Ist die Resp onsevariable bin
ar ko diert und
Bernoulli verteilt, d. h.
y
B
(1
;
), so wird im Logit-Mo dell f
ur den Er-
wartungswert
E
(
y
j
x
)=
P
(
y
=1
j
x
)=
als Resp onsefunktion die logistische
Verteilungsfunktion gew
ahlt. Somit folgt mit dem linearen Pr
adiktor
=
z
0
der Mo dellansatz
=
h
(
)=
exp(
)
1+exp(
)
.
Im ungruppierten Fall liegen N Beobachtungen (
y
i
;x
i
) f
ur
i
=1
;:::;N
vor.
Sind die Daten gruppiert, d. h. zu einem Kovariablenvektor liegen mehrere
unabh
angige Beobachtungen vor, so liegen die Daten in der Form (
y
i
;x
i
;N
i
)
vor. Dab ei ist
y
i
der Resp onse-,
x
i
der Kovariablenvektor und
N
i
die Anzahl
der Beobachtungen in Grupp e
i
. Insgesamt liegen in
N
Grupp en
N
:=
P
N
i
=1
N
i
Beobachtungen vor.
Damit der Mo dellansatz f
ur gruppierte und f
ur ungruppierte Daten mit der
gleichen Notation verwendet werden kann, setzen wir voraus, da die Re-
sp onsevariable
y
i
die relativeH
augkeit in Grupp e
i
b ezeichnet und denie-
ren
w
i
= 1 f
ur ungruppierte Daten und
w
i
=
N
i
f
ur gruppierte Daten. F
ur
die Varianz gilt
var
(
y
i
j
x
i
)=
i
(1
i
)
w
i
, mit
i
=
E
(
y
i
j
x
i
).
Multikategoriales und kumulatives Logit-Mo dell
Diese Mo delle sind
multivariate Erweiterungen des ob en b eschrieb enen Binomialen Logit-Mo dells.
Wir gehen zun
achst von ungruppierten Daten aus, d. h. es ist
y
M
(1
;
)
mit
=(
1
;:::;
q
). Die Auspr
agungen der zugeh
origen Resp onsevariablen
y
i
=(
y
i
1
;:::;y
iq
) liegen in
K
:=
q
+ 1 Kategorien vor und es gilt:
y
ir
=
(
1 falls
Y
i
=
r
0 sonst.
F
ur gruppierte Daten setzen wir voraus, da
y
i
f
ur die einzelnen Kategorien
die Grupp endurchschnitte enth
alt und denieren die Gewichte
w
i
f
ur jede
8
Beobachtungsgrupp e
i
analog zum Binomialen Mo dell.
Im multikategorialen Logit-Mo dell wird der linearen Pr
adiktor mit
i
=
0
B
B
@
i
1
.
.
.
iq
1
C
C
A
=
Z
i
=
2
6
6
6
6
6
6
6
6
4
1
x
0
i
0 0
0
0 0 1
x
0
i
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
0
0
0 1
x
0
i
3
7
7
7
7
7
7
7
7
5
0
B
B
B
B
B
B
B
B
B
B
B
B
@
10
1
20
2
.
.
.
q
0
q
1
C
C
C
C
C
C
C
C
C
C
C
C
A
und die Resp onsefunktion
h
=(
h
1
;:::;h
q
)mit
h
r
(
i
1
;:::;
iq
)=
exp(
ir
)
1+
P
q
s
=1
exp(
is
)
f
ur
r
=1
;:::;q
gew
ahlt.
Sind die Kategorien von
Y
2 f
1
;:::;q
+1 =:
K
g
ordinal, so kann die Ord-
nung der Resp onsekategorien ausgenutzt werden.
Wir b eschr
anken uns auf die Denition f
ur das kumulative Logit-Mo dell.
Weitere Mo delle, welche die Ordnung der Resp onsekategorien ber
ucksichti-
gen, sind b ei Agresti (1990), McCullagh (1980) und Fahrmeir & Tutz (1994)
beschrieb en.
Motiviert wird das kumulative Logit-Mo dell, das auch als Schwellenwert-
mo dell b ezeichnet wird, durch die Existenz einer latenten stetigen Variablen
V
=
x
0
+
, so da mit den Schwellen
1
=
0
<
1
<
<
q
<
K
=
1
Y
=
r
,
r
1
<V
r
f
ur r=1,
:::
,K
9
gilt.
Mit dem Mo dellansatz
P
(
y
r
j
x
) =
F
(
r
+
x
0
), wob ei
F
die logistische
Verteilungsfunktion, d. h.
F
(
u
)=
exp(
u
)
1+exp(
u
)
ist, sowie dem linearen Pr
adiktor
=
0
B
B
@
1
.
.
.
q
1
C
C
A
=
Z
=
2
6
6
6
6
6
6
6
4
1 0
::: :::
0
x
0
0 1 0
:::
0
x
0
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
0
:::
0 1 0
x
0
0
::: :::
0 1
x
0
3
7
7
7
7
7
7
7
5
0
B
B
B
B
B
@
1
.
.
.
q
1
C
C
C
C
C
A
folgt
P
(
Y
=1
j
x
)=
F
(
1
+
x
0
)
P
(
Y
=
r
j
x
)=
F
(
r
+
x
0
)
F
(
r
1
+
x
0
) f
ur
r
=2
;:::;q :
Da f
ur die Schwellen
1
q
erf
ullt sein mu, wird f
ur die Sch
atzung
der Parameter meist eine Parametrisierung der Schwellen mit
1
:=
1
und
r
:= ln(
r
r
1
) f
ur
r
= 2
;:::;q
durchgef
uhrt. Das heit es wird der
Mo dellansatz
=
Z
mit
Z
=
0
B
B
B
B
B
@
1 0
:::
0
x
0
0 1 0
.
.
.
.
.
.
0 1 0
1
C
C
C
C
C
A
und
0
=(
1
;:::;
q
;
)
:
verwendet.
Die Reparametrisierung der Schwellen erfolgt mit
1
=
1
und
r
=
1
+
P
r
s
=2
exp(
s
) f
ur
r
=2
;:::;q
.
10
In b eiden Mo dellen gilt f
ur die Kovarianz
cov
(
y
i
)=
1
w
i
(
diag
(
i
)
i
0
i
), wo-
bei
w
i
und
i
analog zum binomialen Logit-Mo dell deniert sind.
1.1.2 Sch
atzung der Parameter
Der unbekannte Parameter
= (
1
;:::;
p
) wird aus dem Mo dellansatz
E
(
y
i
j
x
i
)=
i
=
h
(
i
) mit
i
=
Z
i
b estimmt. Die Sch
atzung des Parame-
ters kann durch Maximierung der Log-Likeliho o d
l
(
) durchgef
uhrt werden.
Dab ei wird die Nullstelle der Scorefunktion
s
(
) =
@l
@
b erechnet, wob ei zu
beachten ist, da die b eobachtete Fishermatrix
F
obs
(
)=
@
2
l
(
)
@ @
0
p ositiv de-
nit ist.
Da vorausgesetzt wurde, da die
y
i
unabh
angig, identischverteilt sind und die
Verteilung aus der Exp onentialfamilie stammt, gilt (bis auf einen konstanten
Term):
l
(
)=
N
X
i
=1
l
i
(
)=
N
X
i
=1
ln
f
(
y
i
j
i
;;w
i
)=
N
X
i
=1
y
i
i
b
(
i
)
w
i
und f
ur die Scorefunktion folgt mit
1:
s
(
)=
@l
@
=
N
X
i
=1
s
i
(
)=
N
X
i
=1
Z
0
i
D
i
(
)
1
i
(
)[
y
i
i
(
)]
mit
D
i
(
)=
@h
(
i
)
@
und
i
(
)=
cov
(
y
i
).
Die numerische Berechnung des Maximum Likeliho o d Sch
atzers kann unter
Verwendung des Fisher-Scoring Algorithmus b estimmt werden. Ist
F
(
) =
cov
(
s
(
)) =
P
i
F
i
(
)die Fisherinformationsmatrix, so gilt mit
F
(
) =
E
(
F
obs
(
)):
F
i
(
)=
Z
0
i
D
i
(
)
1
i
(
)
D
i
(
)
0
Z
i
:
11
Nach der Initialisierung des Startwertes
^
h
0
i
wird die Iteration
^
h
k
+1
i
=
^
h
k
i
+
F
1
^
h
k
i
s
^
h
k
i
f
ur
k
=0
;
1
;
2
;:::
so lange durchgef
uhrt bis das Abbruchkriterium
jj
^
h
k
+1
i
^
h
k
i
jj
jj
^
h
k
i
jj
<
f
ur ein vorgegeb enes
erf
ullt ist.
12
Kapitel 2
Mo dellerweiterungen im
Generalisierten Linearen
Mo dell
Im vorherigen Kapitel wurde das grundlegende Generalisierte Lineare Mo-
dell vorgestellt. Bei den Mo dellannahmen wurde von b edingt unabh
angigen
Beobachtungen ausgegangen.
Sind die zugrundeliegenden Daten in einer Studie erhob en worden, bei der
einzelne Individuen mehrfach untersucht wurden, so kann keinesfalls davon
ausgegangen werden, da die Daten (
y
i
;x
i
) f
ur
i
= 1
;:::;N
bedingt un-
abh
angig sind. Ab er auch wenn die einzelnen Untersuchungseinheiten zu
Grupp en zusammengefat sind, mu diese Grupp enstruktur durchdas Mo-
dell ber
ucksichtigt werden. Dies ist b eispielsweise dann der Fall, wenn sich
die Untersuchungseinheiten in gemeinsame Familien o der sonstige soziale, re-
gionale o der administrative Grupp en gliedern lassen.
Wir werden diese Grupp enstruktur durch die Mo dellierung von zuf
alligen
Eekten b er
ucksichtigen. Durch diesen Ansatz b er
ucksichtigen wir nichtnur
13
die Grupp enstruktur, sondern es k
onnen zus
atzlich die Ein
usse von nicht
erhob enen Kovariablen untersucht werden.
Diese grupp en- und/o der individuensp ezischen Eekte k
onnen auftreten,
wenn die Erhebung aller Daten aus technischen bzw. wirtschaftliche Gr
unden
nichtm
oglichwar o der wichtige Variablen aufgrund von fehlerhaften Versuchs-
pl
anen
ub ersehen wurden.
Wurden die Daten einer Studie in einem Zeitverlauf gewonnen, so mu zu-
s
atzlich diese Zeitabh
angigkeit b er
ucksichtigt werden. Wir b etrachten hierzu
Mo delle, die eine Variation der Mo dellparameter
ub er die Zeit erm
oglichen.
In diesem Kapitel werden wir daher im folgenden Mo dellans
atze, die entweder
zeitvariierende Ko eÆzienten oder zuf
allige Eekte ber
ucksichtigen, vorstel-
len.
Die Mo dellans
atze f
ur variierende Ko eÆzienten werden in der Literatur in
sehr untschiedlichen Ans
atzen betrachtet. F
ur die Betrachtung der Variati-
on der Mo dellparametern im Generalisierten Linearen Mo dell als Gl
attungs-
funktion von anderen Variablen werden in Hastie & Tibshirani (1993) ver-
schiedene Mo delle vorgestellt. Aus dieser Mo dellklasse werden in Tutz &
Kauermann (1997) lokal gewichtete Sch
atzungen b etrachtet.
In Friedman (1991) und Stone, Hansen, Ko op erb erg & Truong (1997) werden
Sch
atzverfahren mit Gl
attung durch Splines vorgestellt, wob ei insb esondere
die Wahl der Knoten f
ur die Splines b er
ucksichtigt wird. Allgemeine Ans
atze
mit Spline-Funktionen werden in Eubank (1988) sowie in Eilers & Marx
(1996) f
ur B-splines untersucht.
In Silverman (1984) werden sowohl Gl
attungsverfahren mit Spline-Funktionen
als auch Kerngl
attungen mit exibler Bandbreite vorgestellt. Hastie & Loa-
der (1993) schlagen eine lokale Regression mit Kern-Gl
attung vor. Parame-
trische und semi-parametrische Gl
attungsans
atze in Verweildauermo dellen
14
werden in Efron (1988) b etrachtet. Die Zeitabh
angigkeit bei Mewieder-
holungen, die b eispielsweise bei der Analyse von klinischen Studien auftritt,
untersucht Davis (1991) durch semi-parametrische und nichtparametrische
Ans
atze.
Die in dieser Arb eit verwendeten Verfahren f
ur die Sch
atzung der Parameter
mit normalverteilten zuf
alligen Eekten werden in Fahrmeir & Tutz (1994),
Hennevogl (1991) und Tutz & Hennevogl (1996) vorgestellt.
Generalisierte Lineare Mischmo delle mit unterschiedlichen Sch
atzmetho den
werden in McCullo ch (2000) vorgestellt. In Geyer & Thompson (1992) wer-
den Approximierungstechniken b ehandelt, bei denen die Likeliho o dfunktion
direkt maximiert wird. In McCullo ch (1994), McCullo ch (1997) und Bo oth
&Hob ert (1999) werden Verfahren vorgestellt, die einen EM-Algorithmus
mit Gibbs-Sampling bzw. Metrop olis-Hastings Algorithmus verwenden. Ein
anderer Ansatz wird in Rupp ert, Reish & Carroll (1984) und Quintana, Lui
& Pino (1999) gew
ahlt, es wird ein sto chastischer Approximationsalgorith-
mus verwendet und die Likelho o d-Gleichung wird als Regressions Problem
betrachtet.
Da in Studien meist sowohl grupp ensp ezische Eekte als auch Ein
usse
ub er die Zeit b er
ucksichtigt werden m
ussen, werden wir im n
achsten Kapitel
b eide Mo dellans
atze in ein gemeinsames Mo dell aufnehmen. Dieser Ansatz
wird in der Literatur, bis auf Tutz (1999), nicht weiter b etrachtet.
15
2.1 Mo dellans
atze mit zuf
alligen Eekten
Liegen die Beobachtungen von verschiedenen Grupp en bzw. Kategorien, die
als Cluster b ezeichnet werden, vor, kann der Einfu zwischen dem Resp onse
und den erkl
arenden Variablen von Cluster zu Cluster variieren. Deshalb
mu die Grupp enstruktur der Daten im Mo dellansatz ber
ucksichtigt wer-
den. Insb esondere bei Mewiederholungen an mehreren Untersuchungsein-
heiten in einem Zeitverlauf, die im n
achsten Abschnitt b etrachtet werden, ist
mit dieser sogenannten \Extra Variation innerhalb der Daten", die auch als
unb eobachtete Populationsheterogenit
at b ezeichnet wird, zu rechnen.
Da die Beobachtungen nicht in unabh
angigen elementaren Beobachtungsein-
heiten vorliegen, ist die Mo dellannahme im Generalisierten Linearen Mo dell,
welche von gleicher Varianzstruktur f
ur alle Beobachtungen ausgeht, ver-
letzt, und es k
onnen Mo delle, die eine endliche Mischungsdichte zulassen,
angewendet werden. Diese Mo delle werden zum Beispiel in Everitt & Hand
(1981), Laird & Louis (1982), Little & Rubin (1987), McLachlan & Basford
(1988), Meng & van Dyk (1997b) und Titterington, Smith & Makov (1985)
vorgestellt.
Ist die Anzahl der Cluster sehr klein, so ist es m
oglich mit Hilfe von ANOVA
(analysis of variances) und ANCOVA (analysis of covariances) Mo dellen die
Variation zwischen den Clustern zu untersuchen. Man erh
alt Aussagen
ub er
die Variation in den konkret vorliegenden Clustern.
Wenn die Auswahl der b eobachteten Cluster als zuf
allig angesehen wird, ist
man prim
ar nicht an einer Aussage
ub er die konkret vorliegenden Cluster,
sondern
ub er die zugrunde liegende Population interessiert. Es k
onnen Mo-
delle mit zuf
alligen Eekten, die b eispielsweise in Fahrmeir & Tutz (1994),
Jansen (1993), Longford (1993), Meng & van Dyk (1997a), Meng (1997) und
16
Hennevogl (1991) vorgestellt werden, angewendet werden.
Im Mo dellansatz mit zuf
alligen Eekten wird davon ausgegangen, da die
Unterschiede, die in der Zusammenhangsanalyse in den Clustern auftreten,
durch Variablen erkl
art werden k
onnten, die jedo ch nicht b eobachtet wur-
den. Diese unb eobachteten Variablen k
onnen b eispielsweise auftreten, wenn
ihre Erhebung aus technischen bzw. wirtschaftlichen Gr
unden nichtm
oglich
war o der diese relevante Variablen aufgrund von fehlerhaften Versuchspl
anen
ub ersehen wurden.
Im Generalisierten Linearen Mo dell mit zuf
alligen Eekten werden die nicht
b eobachteten Variablen durchzus
atzliche Parameter im linearen Pr
adiktor
ber
ucksichtigt, die als Realisation einer Zufallsvariable aufgefat werden.
In der Denition f
ur Generalisierte Lineare Mo delle mit zuf
alligen Eekten
setzen wir voraus, da f
ur
i
= 1
;:::;N
Beobachtungseinheiten jeweils
n
i
Mewiederholungen vorliegen. Somit liegen insgesamt
N
=
P
N
i
=1
n
i
Beob-
achtungen vor, die sich jedo ch auf
i
= 1
;:::;N
Grupp en, bzw. allgemeiner
ausgedr
uckt auf Beobachtungseinheiten verteilen. Wir setzten damit die Da-
tenstruktur (
y
it
;x
it
)mit
i
=1
;:::;N
und
t
2f
t
i
1
;:::;t
in
i
g
voraus.
Im anschlieenden Kapitel werden wir den Ansatz der zeitvariierenden Ko-
eÆzienten mit der Mo dellierung von zuf
alligen Eekten verbinden. Dort
werden wir voraussetzen, da die Beobachtungen (
y
it
;x
it
)der
i
-ten Beobach-
tungseinheit zum Zeitpunkt
t
erhob en wurden.
Im folgenden b etrachten wir Mo delle b ei denen f
ur jede Beobachtungseinheit
i
ein zuf
alliger Eekt
b
i
in den linearen Pr
adiktor aufgenommen wird. Das
heit, statt dem linearen Pr
adiktor
it
=
Z
it
wird der lineare Pr
adiktor
it
=
Z
it
+
W
it
b
i
b etrachtet. Die Designmatrizen
Z
it
und
W
it
werden durch
geeignete Abbildungen aus den Kovariablen gebildet.
17
Eine
ubliche Annahme f
ur die zuf
alligen Eekte
b
1
;:::;b
N
ist, da diese
unabh
angig normalverteilt sind. Beispiele f
ur diesen Mo dellansatz werden
wir in den Abschnitten 2.1.1 und 2.1.2 b etrachten. Bei dem in Abschnitt
2.1.3 vorgestellten Verfahren der nichtparametrischen Maximum Likeliho o d
wird keine parametrische Voraussetzung f
ur die zuf
alligen Eekte getroen.
Mo delldenition
Die Mo delldenition ist zweistug und teilt sich in das Beobachtungs- und
das Wahrscheinlichkeitsmo dell auf.
Das Beobachtungsmo dell setzt voraus, da die Resp onsewerte
y
i
1
;:::;y
in
i
f
ur
i
=1
;:::;N
b edingt unabh
angig sind, gegeb en
x
i
1
;:::;x
in
i
und
b
i
sowie die
Dichte der
y
it
zu der Exp onentialfamilie geh
ort. Der Erwartungswert
it
ist
ub er die Linkfunktion
h
mit dem linearen Pr
adiktor durch
it
=
E
(
y
it
j
b
i
;x
it
)=
h
(
it
)mit
it
=
Z
it
+
W
it
b
i
verbunden.
Das Wahrscheinlichkeitsmo dell b estimmtdieVerteilung der zuf
alligen Eekte
b
1
;:::;b
N
.In den Abschnitten 2.1.1 und 2.1.2 setzten wir voraus, da diese
f
ur
i
=1
;:::;N
unabh
angig normalverteilt sind mit
b
i
iid
N
(0
;Q
)
:
Die zu sch
atzenden Mo dellparameter in diesem Mo dellansatz sind der Pa-
rametervektor
,der den globalen Einu der Regressoren beschreibt, die
Kovarianzmatrix
Q
und die zuf
alligen Eekte
b
i
, welche die clustersp ezi-
schen Abweichungen vom globalen Einu
beschreib en.
18
F
ur die Sch
atzung der Parameter mit normalverteilten zuf
alligen Eekten
geb en wir im folgenden zwei unterschiedliche Verfahren an. Das erste Ver-
fahren ist der EM-Algorithmus mit Gau-Hermite Integrationstechnik und
das zweite Verfahren der EM-Typ e. Beide Verfahren werden in Fahrmeir &
Tutz (1994), Hennevogl (1991) und Tutz & Hennevogl (1996) vorgestellt.
2.1.1 EM-Algorithmus mit Gau-Hermite Integrations-
technik
Die Sch
atzung der Parameter in diesem Algorithmus erfolgt in zwei Schritten.
Zuerst werden die festen Parameter
und
Q
gesch
atzt und anschlieend die
zuf
alligen Eekte
b
i
.
Sch
atzung der festen Parameter
Um die Notation
ub ersichtlicher zu gestalten, b eschreib en wir kurz das ver-
wendete Verfahren f
ur einen eindimensionalen zuf
alligen Eekt. Im Anhang
im Abschnitt 6.1.2 gehen wir auf die Details des Verfahrens ein und b etrach-
ten dab ei den allgemeinen Fall b ei dem die zuf
alligen Eekte
b
i
mehrdimen-
sional sein k
onnen. Im folgenden eindimensionalen Fall schreib en wir statt
der Kovarianzmatrix
Q
einfach
2
.
Der lineare Pr
adiktor im univariaten Generalisierten Linearen Mo dell mit
einem eindimensionalen zuf
alligen Eekt
it
=
x
0
it
+
b
i
mit
b
i
iid
N
(0
;
2
)
l
at sich mit
b
i
=
a
i
umparametrisieren zu
it
=
x
0
it
+
a
i
mit
a
i
iid
N
(0
;
1)
:
19
Mit der
ublichen Voraussetzung der b edingten Unabh
angigkeit zwischen den
einzelnen Beobachtungseinheiten bzw. innerhalb der Mewiederholungen zu
einer Beobachtungseinheit gilt f
ur die Likeliho o d
L
(
;
)=
N
Y
i
=1
Z
f
(
y
i
j
a
i
;
;
)
p
(
a
i
)
da
i
mit
f
(
y
i
j
a
i
;
;
)=
n
i
Y
t
=1
f
(
y
it
j
a
i
;
;
)
:
Wob ei
f
(
) die b edingte Dichte, die aus einer Exp onentialfamilie stammtund
p
(
) die Dichte der Standardnormalverteilung ist.
Diese Mischungsdichte kann nur f
ur einige Sp ezialf
alle analytischgel
ost wer-
den, so da f
ur die L
osung des Integrals numerische Verfahren angewendet
werden.
Bei der Gau-Hermite Integrationstechnik wird das Integral durch eine endli-
che Summe
ub er M Quadratur-St
utzstellen
d
j
,(
j
=1
;:::;M
)mitbekannten
Quadratur-Gewichten
v
j
,
j
=1
;:::;M
approximiert, so da
L
(
;
)
N
Y
i
=1
M
X
j
=1
v
j
f
(
y
i
j
d
j
;
;
)
folgt. Die Gewichte
v
j
und St
utzstellen
d
j
f
ur die Gau-Hermite Integration
liegen in Tab ellen (Abramowitz & Stegun, 1972) vor.
Das b edeutet, da die Likeliho o d durch eine endliche Mischung von Dichten
aus einer Exp onentialfamilie mit b ekannten Mischungsgewichten
v
j
und be-
kannten Mapunkten
d
j
approximiert wird. Ist
l
(
;
)die Log-Likelihood,
d. h.
l
(
;
)=
N
X
i
=1
ln
M
X
j
=1
v
j
f
ij
mit
f
ij
=
f
(
y
i
j
d
j
;
;
)
20
so folgt (Hinde, 1982; Aitkin, 1996)
@l
@
=
N
X
i
=1
M
X
j
=1
v
j
f
ij
@
ln
f
ij
@
M
P
l
=1
v
l
f
il
=
N
X
i
=1
M
X
j
=1
c
ij
s
ij
(
)
@l
@
=
N
X
i
=1
M
X
j
=1
v
j
f
ij
@
ln
f
ij
@
M
P
l
=1
v
l
f
il
=
N
X
i
=1
M
X
j
=1
c
ij
s
ij
(
)
:
Dab ei ist
c
ij
die p osteriori Wahrscheinlichkeit, da die Beobachtungseinheit
i
aus der Mischungskomp onente
j
stammt, mit
c
ij
=
v
j
f
ij
M
P
s
=1
v
s
f
is
:
Die Komp onenten
s
ij
(
)und
s
ij
(
)sind die
-bzw.
-Komp onenten der
Scorefunktion f
ur die
i
-te Beobachtung in der
j
-ten Mischungskomp onente.
Da die St
utzstellen
d
j
in Tab ellen vorliegen, sind diese Gr
oen b ekannt. Die
Parameter
und
sind die zugeh
origen Regressionsko eÆzienten.
Beim Maximieren der Log-Likeliho o d wird die Scorefunktion gleich Null ge-
setzt, und man erh
alt die gleichen gewichteten Summen, die als Gleichun-
gen eines gew
ohnlichen Generalisierten Linearen Mo dells mit dem linearen
Pr
adiktor
itj
=
x
0
it
+
d
j
=[
x
0
it
; d
j
]
!
mit
8
>
>
<
>
>
:
i
=1
;:::;N
t
=1
;:::;n
i
j
=1
;:::;M
entstehen, wenn die Beobachtungen mit
c
ij
gewichtet werden und f
ur den
zugeh
origen Resp onse
y
itj
y
it
gilt.
21
Das L
osen der Gleichungen f
ur gegeb ene Gewichte
c
ij
und die Berechnung der
Gewichte mit der letzten Parametersch
atzung f
ur
0
=(
0
;
)f
uhrt Fahrmeir
&Tutz (1994) zu folgendem EM-Algorithmus:
Die Anzahl
M
der St
utzstellen f
ur die Gau-Hermite Integration wird
schrittweise erh
oht.
Die Daten des Regressionsmo dells werden mit der Anzahl der St
utz-
stellen vervielfacht, indem jede Zeile
y
it
=
x
0
it
des Mo dellansatzes
durch
0
B
B
@
y
it
.
.
.
y
it
1
C
C
A
=
2
6
6
4
x
0
it
d
1
.
.
.
.
.
.
x
0
it
d
M
3
7
7
5
0
B
B
@
1
C
C
A
mit
=
!
ersetzt wird.
Der Startwert des Sch
atzers f
ur die Iteration ist ^
h
0
i
= (
^
h
0
i
;
^
h
0
i
)
0
,
wob ei
^
h
0
i
der Sch
atzer des Generalisierten Linearen Mo dells ohne
zuf
allige Eekte mit dem linearen Pr
adiktor
it
=
x
0
it
ist. F
ur ^
h
0
i
wird eine Konstante b eispielsweise ^
h
0
i
=2 gew
ahlt.
Jeder Iterationsschritt
k
=0
;
1
;
2
;:::
b esteht aus den folgenden Schrit-
ten:
{
Berechnung der Gewichte mit:
c
ij
(
h
k
i
)=
v
j
f
(
y
i
j
d
j
;
h
k
i
)
M
P
s
=1
v
s
f
(
y
i
j
d
s
;
h
k
i
)
{
Sch
atzung von ^
h
k
+1
i
durch das mit
c
ij
(
h
k
i
) gewichtete Regressi-
onsmo dell.
22
Ist das Abbruchkriterium:
jj
^
h
k
+1
i
^
h
k
i
jj
jj
^
h
k
i
jj
<
f
ur ein vorgegeb enes
erf
ullt, so wird die Iteration b eendet.
Durch die Verwendung des EM-Algorithmus ist sichergestellt, da sich in
jedem Iterationsschritt der Wert der Likeliho o d erh
oht (Wu, 1983). Es ist
jedo ch nicht gew
ahrleistet, da das globale Maximum gefunden wird.
Der Startwert
h
0
i
= 0 sollte deshalb vermieden werden. Denn bei
0
=
(
0
;
0), d. h. b eim Mo dellansatz ohne Ber
ucksichtigung von zuf
alligen Eek-
ten, liegt meist ein lokales Maximum vor. In den durchgef
uhrten Simulatio-
nen hat sich gezeigt, da
h
0
i
=2 eine gute Wahl f
ur den Startwert ist.
Um die optimale Anzahl der gew
ahlten St
utzstellenanzahl
M
bei der Gau-
Hermite Integration zu nden, mu der Algorithmus f
ur
M
= 2
;
4
;
6
;:::
solange durchgef
uhrt werden, bis sichder Wert der Likeliho o d nicht mehr
erh
oht.
Sch
atzung der zuf
alligen Eekte
Nachdem die festen Parameter
und
2
,bzw. allgemein die Kovarianzma-
trix
Q
, gesch
atzt wurden, kann die Bayes-Sch
atzung der zuf
alligen Eekte
b
i
erfolgen.
Der unbekannte Parameter
b
i
wurde als Zufallsgr
oe mit
b
i
N
(0
;Q
) vor-
ausgesetzt, so da die Priori-Dichte mit
p
(
b
i
)=
1
(2
)
q=
2
j
Q
j
1
=
2
exp
1
2
b
0
i
Q
1
b
i
als bekannt angenommen wird. Wir sch
atzen
b
i
durch Maximierung der
Posteriori Dichte
p
(
b
i
j
Y
), d. h. wir maximieren
l
p
(
b
i
j
Y
):=ln
p
(
b
i
j
Y
)mit:
23
l
p
(
b
i
j
Y
)=ln
L
(
b
i
j
Y
)+ln
p
(
b
i
)
=ln
L
(
b
i
j
Y
)
1
2
b
0
i
Q
1
b
i
;
wob ei
L
(
b
i
j
Y
) die Likeliho o dfunktion, abh
angig von
b
i
,bei gegeb enden Da-
ten
Y
b ezeichnet und konstante Terme vernachl
at wurden.
Die Bestimmung von
b
i
kann durchAnwendung des Fisher-Scoring Algorith-
mus durchgef
uhrt werden.
Mit dem Startwert
b
h
0
i
i
= 0 wird die Integration
b
h
k
+1
i
i
=
b
h
k
i
i
+
F
1
(
b
h
k
i
i
)
s
(
b
h
k
i
i
)
f
ur
k
=1
;
2
;:::
durchgef
uhrt bis ein Abbruchkriterium erf
ullt ist, wob ei
F
(
)
die Fishermatrix und
s
(
) die Scorefunktion b ezeichnen.
Als Anwendungsb eispiel b etrachten wir im folgenden das Bin
are Logit-Mo dell
b ei dem die nat
urliche Linkfunktion verwendet wird. Liegt ein eindimensio-
naler zuf
alliger Eekt vor und ist
W
it
1, d.h.
l
p
(
b
i
j
Y
)=ln
L
(
b
i
j
Y
)
1
2
Q
1
b
2
i
so folgt f
ur Scorefunktion
s
(
b
i
)=
n
i
X
t
=1
n
it
[
y
it
it
(
b
i
)]
Q
1
b
i
sowie f
ur die Fishermatrix
F
(
b
i
)=
n
i
X
t
=1
n
it
it
(
b
i
)(1
it
(
b
i
)) +
Q
1
mit
it
(
b
i
) =
h
(
Z
it
+
b
i
), wob ei f
ur
und
Q
die entsprechenden Sch
atzer
verwendet werden.
24
2.1.2 EM-Typ e Algorithmus
Beim EM-Typ e Algorithmus wird eine simultane Sch
atzung der Parame-
ter
,
Q
und
b
i
,
i
= 1
;::: N
mit Hilfe des EM-Algorithmus durchgef
uhrt.
Das grundlegende Prinzip des EM-Typ e Algorithmus wird auchbeim EM-
Algorithmus mit Gau-Hermite Integrationstechnik verwendet und im An-
hang in Abschnitt 6.1.2 vorgestellt.
Kurz zusammengefat kann der EM-Typ e Algorithmus wie folgt b eschrieb en
werden. Die Parameter (
0
;b
0
1
;:::;b
0
N
) =
Æ
0
werden als unbekannte Daten
betrachtet und
Q
ist der zu sch
atzende Parameter.
Mit diesem Ansatz k
onnen im
h
k
+1
i
-ten Iterationsschritt im M-Schritt die
Sch
atzer f
ur
Æ
0
h
k
+1
i
=(
0
h
k
+1
i
;b
0
h
k
+1
i
1
;:::;b
0
h
k
+1
i
N
), gegeb en
Q
h
k
i
b estimmtwerden. Im E-Schritt wird dann
Q
mit dem Ansatz
M
Q
j
Q
h
k
i
=
E
ln
f
(
Y; Æ;Q
)
j
Y; Q
h
k
i
gesch
atzt.
Auf die Einzelheiten des EM-Typ e Algorithmus wird, unterteilt in den E-
bzw. M-Schritt, im folgenden eingegangen.
M-Schritt
Mit gegeb enem
Q
h
k
i
gilt f
ur die p osteriori Log-Likeliho o d
l
(
Æ
)
l
(
Æ
) =
N
X
i
=1
ln
f
(
y
i
j
b
i
;
)
1
2
N
X
i
=1
b
0
i
Q
h
k
i
1
b
i
=
N
X
i
=1
n
i
X
t
=1
ln
f
(
y
it
j
b
i
;
)
1
2
N
X
i
=1
b
0
i
Q
h
k
i
1
b
i
;
25
so da f
ur die Scorefunktion
@l
(
Æ
)
@
=
N
X
i
=1
n
i
X
t
=1
Z
0
it
D
it
1
it
(
y
it
it
)
@l
(
Æ
)
@b
i
=
T
X
t
=1
W
0
it
D
it
1
it
(
y
it
it
)
Q
h
k
i
1
b
i
mit
D
it
=
D
it
(
; b
i
) =
@h
(
it
)
,
1
it
=
1
it
(
; b
i
) = (
cov
(
y
it
j
; b
i
))
1
und
it
=
it
(
; b
i
)=
h
(
it
) folgt. F
ur die Fishermatrix gilt
F
(
Æ
)=
2
6
6
6
6
6
4
F
F
b
1
F
b
N
F
b
1
F
b
1
b
1
O
.
.
.
.
.
.
F
b
N
O
F
b
n
b
n
3
7
7
7
7
7
5
mit
F
=
N
X
i
=1
n
i
X
t
=1
Z
0
it
D
it
1
it
D
0
it
Z
it
F
b
i
=
n
i
X
t
=1
Z
0
it
D
it
1
it
D
0
it
W
it
F
b
i
b
i
=
n
i
X
t
=1
W
0
it
D
it
1
it
D
0
it
+
Q
h
k
i
1
F
b
i
=
F
0
b
i
:
Die Berechnung von
Æ
h
k
+1
i
kann durchden Fisher-Scoring Algorithmus mit
Æ
h
k
+1
i
=
Æ
h
k
i
+
F
1
Æ
h
k
i
s
Æ
h
k
i
26
durchgef
uhrt werden.
Im E-Schritt werden die Sch
atzer
^
b
i
sowie ihre p osteriori Kovarianzen
V
ii
f
ur
i
=1
;:::;n
ben
otigt. Letztere ergeb en sich aus
F
1
(
Æ
)=
2
6
6
6
6
6
4
V
V
1
V
n
V
1
V
11
V
1
n
.
.
.
.
.
.
V
n
V
n
1
V
nn
3
7
7
7
7
7
5
mit
V
=(
F
n
X
i
=1
F
i
F
1
ii
F
i
)
1
V
ii
=
F
1
ii
+
F
1
ii
F
i
V
F
i
F
1
ii
V
i
=
V
0
i
=
V
F
i
F
1
ii
V
ij
=
V
0
ji
=
F
1
ii
F
i
V
F
j
F
1
jj
;i
6
=
j
Im Anhang in Abschnitt 6.2.1 wird gezeigt, da die Blo ckdiagonalgestalt der
Fishermatrix ausgenutzt werden kann, so da f
ur die Berechnung von
Æ
h
k
+1
i
nicht die gesammte Fishermatrix
F
(
Æ
h
k
i
)invertiert werden mu.
E-Schritt
Da f
ur die gemeinsame Dichte
f
(
Y; Æ;Q
)=
N
Y
i
=1
f
(
Y; Æ
)
p
(
b
i
;Q
)
27
gilt und nur der letzte Faktor in der Dichte von
Q
abh
angt, l
at sichder
E-Schritt im
h
k
+1
i
-ten Iterationsschritt reduzieren zu:
~
M
Q
j
Q
h
k
i
=
E
ln
p
(
Æ;Q
)
j
Y; Q
h
k
i
Da
b
i
N
(0
;Q
)vorausgesetzt wurde, erh
alt man bis auf konstante Terme
~
M
Q
j
Q
h
k
i
=
N
2
ln
j
Q
j
N
X
i
=1
tr
Q
1
E
b
i
b
0
i
j
Y; Q
h
k
i
und damit kann
Q
h
k
+1
i
b estimmtwerden mit
Q
h
k
+1
i
=
1
N
N
X
i
=1
cov
(
b
i
j
y
i
;Q
h
k
i
)+
E
(
b
i
j
y
i
;Q
h
k
i
)
E
(
b
i
j
y
i
;Q
h
k
i
)
0
:
In Tutz & Hennevogl (1996) wird vorgeschlagen, wegen des groen nume-
rischen Aufwands bei der Berechnung der Erwartungswerte, den p osteriori
Erwartungswert von
b
i
durchden p osteriori Mo de
^
b
h
k
i
i
zu ersetzen und die
p osteriori Kovarianzen durch
^
V
ii
. Damit folgt
Q
h
k
i
=
1
N
N
X
i
=1
^
V
h
k
i
ii
+
^
b
h
k
i
i
^
b
h
k
i
i
0
:
2.1.3 Mo dellansatz der nichtparametrischen
Maximum Likeliho o d
Im vorherigen Abschnitt wurde vorausgesetzt, da die eindimensionalen zu-
f
alligen Eekte unabh
angig normalverteilt mit
b
i
iid
N
(0
;
2
) sind. Das heit,
es wird eine parametrische Form vorausgesetzt. Da die Mischungsdichte nicht
bekannt ist und die Parametersch
atzungen empndlichvon der Sp ezizierung
der parametrischen Form abh
angen kann, bietet ein nichtparametrischer An-
satz, wie die nichtparametrische Maximum Likeliho o d-Sch
atzung Vorteile.
28
Bei diesem Verfahren, das vorgestellt wird in Aitkin (1995), Aitkin (1996),
Aitkin & Aitkin (1996), Aitkin & Francis (1995) und Aitkin & Wilson (1980),
wird die Mischverteilung zusammen mit den Mo dellparametern des Genera-
lisierten Linearen Mo dells gesch
atzt. Das prim
are Ziel ist jedo ch nicht die
Mischverteilung zu sch
atzen, sondern es sollen falsche Schlufolgerungen auf-
grund fehlerhafter Mo dellvoraussetzungen vermieden werden.
Bei der Parametersch
atzung mit der Gau-Hermite Integration sind die St
utz-
stellen
d
j
und die zugeh
origen Gewichte
v
j
f
ur eine gew
ahlte Anzahl von
St
utzstellen bei der gesamten Iteration konstant. Es wird ein gewichtetes
Generalisiertes Lineares Mo dell mit dem linearen Pr
adiktor
itj
=
x
0
it
+
d
j
betrachtet.
An Stelle von
d
j
werden b ei der nichtparametrischen Maximum Likeliho o d
Sch
atzung die St
utzstellen
1
;:::;
M
als Parameter in den linearen Pr
adiktor
mit
itj
=
x
0
it
+
j
aufgenommen, so da
j
der Parameter f
ur den Interzept der
j
-ten Komp o-
nente der Mischverteilung ist.
Die Parameter
1
;:::;
M
k
onnen gesch
atzt werden, indem Indikatoren mit
M
Stufen in das Mo dell aufgenommen werden. Da Identizierungsprob eme
mit dem Interzept von
und einer Komp onente
j
entstehen, wird
M
0
gesetzt. Alternativ ist es auch m
oglich den Interzept aus dem Mo dell zu
entfernen.
Die Gewichte zu den entsprechenden St
utzstellen
p
1
;:::;p
M
mit
P
M
j
=1
p
j
=1
m
ussen in jedem Iterationsschritt b estimmtwerden. Durch Dierenzieren der
29
Log-Likeliho o d
l
(
;
1
;:::;
M
)=
N
X
i
=1
ln
M
X
j
=1
p
j
f
ij
mit
f
ij
=
f
(
y
i
j
;
1
;:::;
M
)
:
nach
p
j
folgt
@l
@p
k
=
N
X
i
=1
f
ij
f
iM
P
M
l
=1
p
l
f
il
=
N
X
i
=1
w
ik
p
k
w
iM
p
M
und mit
@l
@p
k
!
= 0 gilt:
^
p
j
=
N
X
i
=1
w
ij
N
:
Die Sch
atzung der Parameter kann dadurch mit folgendem EM-Algorithmus
durchgef
uhrt werden:
Die Anzahl
M
der St
utzstellen wird fest gew
ahlt und schrittweise erh
oht
bis das Maximum der Likeliho o d erreicht ist.
Die Daten des Regressionsmo dells werden mit der Anzahl der St
utz-
stellen vervielfacht, indem jede Zeile
y
it
=
x
0
it
des Mo dellansatzes
durch
0
B
B
@
y
it
.
.
.
y
it
1
C
C
A
=
2
6
6
6
6
6
6
6
4
x
0
it
1 0
0
.
.
.
.
.
.
.
.
.
.
.
.
x
0
it
0
:::
0 1
x
0
it
0
:::
0 0
3
7
7
7
7
7
7
7
5
0
B
B
B
B
B
B
B
@
z
1
z
2
.
.
.
z
M
1
1
C
C
C
C
C
C
C
A
ersetzt wird. Es werden jeweils
M
Zeilen, ab er nur
M
1neue Para-
meter in das Regressionsmo dell aufgenommen, so da
M
0 erf
ullt
ist und keine Identizierungsprobleme auftreten.
30
Startwert f
ur die Iteration ist die Berechnung der Gewichte
w
h
0
i
ij
,wob ei
h
0
i
der Sch
atzer des Generalisierten Linearen Mo dells ohne zuf
allige
Eekte ist und die St
utzstellen und Gewichte aus den Verfahren mit
Gau-Hermite Integration
ub ernommen werden, so da
z
h
0
i
j
=
d
j
und
p
h
0
i
j
=
v
j
f
ur
j
=1
;:::;M
gilt, d. h. es ist:
w
h
0
i
ij
=
v
j
f
(
y
i
j
h
0
i
;d
j
)
M
P
l
=1
v
l
f
(
y
i
j
h
0
i
;d
l
)
:
Jeder Iterationsschritt
k
=0
;
1
;
2
;:::
b esteht aus den folgenden Schrit-
ten:
{
Sch
atzung von ^
h
k
+1
i
= (
^
h
k
+1
i
;
^
z
h
k
+1
i
1
;:::;
^
z
h
k
+1
i
M
1
)
0
durch das mit
w
ij
(
h
k
+1
i
) gewichtete Regressionsmo dell
{
Up date der p osteriori Wahrscheinlichkeit durch
p
h
k
+1
i
j
=
N
P
i
=1
w
h
k
i
ij
N
{
Berechnung der Gewichte f
ur Beobachtungseinheit
i
und St
utzstelle
j
mit
w
ij
(
h
k
+1
i
)=
p
h
k
i
j
f
(
y
i
j
h
k
i
;z
h
k
i
j
)
M
P
l
=1
p
h
k
i
l
f
(
y
i
j
h
k
i
;z
h
k
i
l
)
:
Die Iteration wird f
ur eine fest gew
ahlte Anzahl von St
utzstellen solan-
ge durchgef
uhrt, bis ein lokales Maximum f
ur die Likeliho o d erreicht
wurde. Dann wird die Anzahl der St
utzstellen erh
ohtbiskein gr
oerer
Wert f
ur das Maximum der Likeliho o d gefunden wird.
31
Damit die Designmatrix des Regressionsansatzes vollen Rang hat, wurde eine
Parametrisierung des Mo dellansatzes
itj
=
x
0
it
+
j
durch:
~
it
=
0
B
B
@
~
it
1
.
.
.
~
itM
1
C
C
A
=
2
6
6
6
6
6
6
6
4
x
0
it
1 0
0
.
.
.
.
.
.
.
.
.
.
.
.
x
0
it
0
:::
0 1
x
0
it
0
:::
0 0
3
7
7
7
7
7
7
7
5
0
B
B
B
B
B
B
B
@
~
~
1
~
2
.
.
.
~
M
1
1
C
C
C
C
C
C
C
A
mit ~
M
0gew
ahlt.
Unter der Neb enb edingung
M
P
j
=1
j
=0 k
onnen die urspr
unglichen Parameter
=(
0
;
1
:::;
q
)und
=(
1
;:::;
M
) mit:
o
=
~
0
+
1
M
M
X
j
=1
~
j
j
=
~
j
f
ur
j
=1
;:::;q
j
=
~
0
0
+~
j
f
ur
j
=1
;:::;M
reparametrisiert werden.
32
2.2 Mo dellans
atze mit (zeit-)variierenden
Ko eÆzienten
Beim Mo dellansatz mit (zeit-)variierenden Ko eÆzienten b etrachten wir den
zeitabh
angigen Regressionsparameter
(
t
) und setzen voraus, da die Daten
in der Form (
y
i
;x
i
;t
) vorliegen. Das heit, da neb en dem Resp onse- und
Kovariablenvektor auch no ch der Beobachtungszeitpunkt
t
f
ur jede Beob-
achtung bekannt ist. Im allgemeinen gehen wir davon aus, da f
ur die
i
-te
Beobachtungseinheit mehrere Beobachtungszeitpunkte bzw. Mewiederho-
lungen zu verschiedenen Zeitpunkten vorliegen.
Weiter wird angenommen, da die Menge der beobachteten Zeitpunkte dis-
kret und endlich ist, so da
t
2f
1
;:::;T
g
vorausgesetzt werden kann.
Wir verwenden die Daten (
y
i
;x
i
;t
) in der Notation (
y
it
;x
it
)mit
i
=1
;:::;N
und
t
2f
1
;:::;T
g
,wob ei
N
die Anzahl der Beobachtungseinheiten ist.
Liegt ein binomiales o der multinomiales Mo dell mit gruppierten Daten vor, so
wird vorausgesetzt, da die Daten in der Form (
y
it
;x
it
;n
it
)mit
i
=1
;:::;N
und
t
2f
1
;:::;T
g
vorliegen. Dab ei b ezeichnet
N
die Anzahl der Grupp en,
n
it
die Anzahl der Beobachtungen,
x
it
den Kovariablenvektor und
y
it
den
Grupp endurchschnitt jeweils f
ur die Grupp e
i
zum Zeitpunkt
t
. Die Anzahl
der Beobachtungen ist
N
=
P
i
P
t
n
it
.
Im weiteren unterscheiden wir nicht zwischen gruppierten und ungruppierten
Daten. Wir setzen
n
it
= 1 f
ur ungruppierte Daten und sprechen allgemein
von der
i
-ten Beobachtungseinheit.
Es wurde nicht vorausgesetzt, da zu einer Beobachtungseinheit
i
f
ur alle
Zeitpunkte
t
=1
;:::;T
Beobachtungen vorliegen. Trotzdem werden wir die
Summation
P
N
i
=1
P
T
t
=1
verwenden, indem wir die entsprechenden Summan-
den mit 0 denieren, wenn f
ur eine Beobachtungseinheit
i
zum Zeitpunkt
t
33
keine Beobachtung vorliegt.
Um Zeitabh
angigkeit zu mo dellieren, b etrachten wir den RegressionskoeÆ-
zienten
als eine Funktion der Zeit
t
,d. h.
0
=
0
(
t
)=(
1
(
t
)
;:::;
p
(
t
)).
Die in dieser Arb eit verwendeten Ans
atze f
ur
j
(
t
) (
j
2
1
;:::;p
) werden
im folgenden kurz vorgestellt. Dab ei gehen wir immer davon aus, da Zeit-
abh
angigkeit vorliegt, auch wenn mit den Mo dellen allgemeinere Abh
angig-
keitsstrukturen mo delliert werden k
onnen.
2.2.1 Parametrische und semi-parametrische Ans
atze
Polynomiale Regression
Bei der Polynomialen Regression wird die Zeitabh
angigkeit des Regessions-
ko eÆzienten
j
(
t
) durch ein Polynom von
t
mit fest vorgegeb enem Grad
s
mo delliert:
j
(
t
)=
s
X
r
=0
jr
t
r
:
Regressions Spline
Bei diesem Ansatz (Hastie & Tibshirani, 1990; Fahrmeir, Hamerle & Tutz,
1994) wird der obige Mo dellansatz ver
andert indem der Parameter
j
(
t
)
st
uckweise durch Polynome deniert wird. Der Spline ist an den Bruch-
stellen, den sogenannten Knoten stetig dierenzierbar.
Eine m
ogliche Anwendung ist der st
uckweise-kubische Spline Efron (1988).
Dab ei wird
j
(
t
) aus
s
+1 (
s
2
N
) kubischen Funktionen zusammengesetzt, so
da die Funktion auch an den Knoten
1
;:::;
s
,zweimal stetig dierenzierbar
34
ist:
j
(
t
)=
j
0
+
j
1
t
+
j
2
t
2
+
j
3
t
3
+
s
X
r
=1
jr
+3
(
t
r
)
3
+
;
mit (
t
r
)
+
=
max
f
0
;t
r
g
.
Liegen die Beobachtungen teilweise sp
arlich vor, b eispielsweise ist dies in
Verweildauermo dellen gegen Ende des Beobachtungszeitraumes der Fall, so
wird h
aug ein kubisch-linearer Regressions Spline verwendet:
j
(
t
)=
j
0
+
j
1
t
+
j
2
(
t
)
2
+
j
3
(
t
)
3
mit (
t
)
=
min
f
0
;t
g
.
F
ur
t<
hat die Funktion
j
(
t
) einen kubischen und f
ur
t>
einen linea-
ren Verlauf. Die Funktion
j
(
t
)ist f
ur alle t, insb esondere f
ur
t
=
stetig
dierenzierbar.
Neb en diesen Regressions Splines ist es auchm
oglich, mehrere kubische Funk-
tionen und eine lineare Funktion st
uckweise zu denieren, so da die Funktion
j
(
t
) f
ur
t>
aus einer Geraden und f
ur
t <
aus einem st
uckweise-
kubischen Spline b esteht, es ist:
j
(
t
)=
j
0
+
j
1
t
+
j
2
(
t
)
2
+
j
3
(
t
)
3
+
s
X
r
=1
jr
+3
(
t
r
)
3
:
Mit den Knoten
r
<
r
1
<:::;
1
<
ist
j
(
t
) f
ur
t
einmal und f
ur
t<
zweimal stetig dierenzierbar.
Die in diesem Abschnitt vorgeschlagenen Mo delle sind sehr einfach zu im-
plementieren. Nach Anpassung der Designmatrix k
onnen Standardverfahren
zur Sch
atzung der Ko eÆzienten
jr
verwendet werden um
j
(
t
)zu b estim-
men.
35
Um die Interpretation der variierenden Ko eÆzienten zu erleichtern und um
numerische Sch
atzprobleme zu vermeiden, mu darauf geachtet werden, da
der Grad des Polynoms bzw. die Anzahl der Knoten der vorhanden Daten-
struktur gerecht wird. Das heit die Anzahl der Spalten in der Designmatrix
sollte m
oglichst klein sein, um die Struktur der zeitlichen Abh
angigkeit zu
erkennen und um eine
Ub erparametrisierung zu vermeiden
2.2.2 Lokaler Likeliho o dansatz
Dieser Mo dellansatz, der auchin Kauermann & Tutz (1995), Kauermann
& Tutz (2000), Tutz & Kauermann (1995) und Tutz & Kauermann (1997)
b ehandelt wird, ist ein sp ezielles variierendes Ko eÆzientenmo dell, wie es in
Hastie & Tibshirani (1993) vorgestellt wird. Die Mo dellparameter variieren
in Abh
angigkeit einer erhob enen Variablen, dem sogenannten Eektmo di-
zierer. In unserem Ansatz entspricht der Eektmo dizierer der Zeit.
Zu jedem Zeitpunkt
t
2
1
;:::;T
wird
(
t
) = (
1
(
t
)
;:::;
p
(
t
))
0
durchein
separates Generalisiertes Lineares Mo dell gesch
atzt. Dab ei gehen die Beob-
achtungen (
y
is
;x
is
)mit
s
2 f
1
;:::;T
g
gewichtet, in Abh
angigkeit zu dem
b eobachteten Zeitpunkt
s
, in die Sch
atzung ein. Die dab ei verwendete Ge-
wichtsfunktion
!
(
t; s
) ist unter Verwendung einer Kernfunktion
K
(
u
) in
Abh
angigkeit des Gl
attungsparameters
durch
!
(
t; s
)=
c
t
K
t
s
deniert.
Als Kernfunktionen werden verschiedene unimo dale, symmetrische Funk-
36
tionen verwendet. Beispielsweise k
onnen der Epanechnikov-Kern
K
(
u
)=
(
3
4
(1
u
2
) f
ur
j
u
j
1
0sonst
o der der Gau-Kern
K
(
u
)=
1
p
2
exp
1
2
u
2
:
angewendet werden.
Mit der Wahl der Normierungskonstanten
c
t
=
K
(0)
1
gilt:
!1
:
!
(
t; s
)=1 f
ur
s
=1
;:::;T
!
0:
!
(
t; s
)=
(
1 f
ur
s
=
t
0 f
ur
s
6
=
t
d. h. im Lokalen Likeliho o dansatz sind auchzwei Sp ezialf
alle enthalten. F
ur
!1
erhalten alle Beobachtungen unabh
angig von Beobachtungszeitpunkt
das Gewicht 1. Die Parametersch
atzung mit dem Lokalen Likeliho o dansatz
entspricht der Sch
atzung mit dem Mo dellansatz b ei dem keine Zeitabh
angig-
keit f
ur
ber
ucksichtigt wird. F
ur
!
0 gehen in die Parametersch
atzung
von
(
t
)nur die Beobachtungen zum Zeitpunkt
t
ein.
Durch die Wahl des Gl
attungsparameters k
onnen somit sehr exib el die un-
terschiedlichen Gewichte festgelegt werden. In Abbildung 2.1 werden die
Gewichte
!
(
t; s
) f
ur
t
=10 und
s
= (1
;
2
;:::;
21) in Abh
angigkeit des
Gl
attungsparameters
unter Verwendung des Gau-Kerns dargestellt.
Liegt ein Verweildauermo dell vor, b ei dem gegen Ende des Beobachtungszeit-
raumes nur no ch wenige Beobachtungen vorliegen, so kann der Gl
attungs-
parameter
in Abh
angigkeit der Anzahl der vorliegenden Beobachtungen
b eispielsweise mit
t
=
1
p
n
t
gew
ahlt werden.
37
5 101520
0.0 0.2 0.4 0.6 0.8 1.0
gamma = 0.1
5101520
0.0 0.2 0.4 0.6 0.8 1.0
gamma = 0.5
5 101520
0.0 0.2 0.4 0.6 0.8 1.0
gamma = 1
5 101520
0.0 0.2 0.4 0.6 0.8 1.0
gamma = 2
5101520
0.0 0.2 0.4 0.6 0.8 1.0
gamma = 3
5 101520
0.0 0.2 0.4 0.6 0.8 1.0
gamma = 4
5 101520
0.0 0.2 0.4 0.6 0.8 1.0
gamma = 5
5101520
0.0 0.2 0.4 0.6 0.8 1.0
gamma = 10
5 101520
0.0 0.2 0.4 0.6 0.8 1.0
gamma = 100
Abbildung 2.1: Gewichte des Gau-Kerns in Abh
angigkeit des
Gl
attungsparameters
Die Wahl des Gl
attungsparameters
kann durch Kreuzvalidierung , die b ei-
spielsweise in Rice & Silverman (1991), Fan, Hall, Martin & Patil (1996) und
Tutz (1989)) b ehandelt wird, erfolgen. Dab ei wird der Wert von
gew
ahlt
b ei dem
CV
(
) mit
CV
(
)=
1
N
N
X
i
=1
T
X
t
=1
l
(
y
it
;
^
(
i;t
)
it
)
minimal ist.
Dab ei ist ^
(
i;t
)
it
der gesch
atzte Resp onse, gegeb en
x
it
, der durch den Sch
atzer
38
d
(
t
)
(
i;t
)
b estimmt wird. Der Sch
atzer
d
(
t
)
(
i;t
)
ist derjenige Sch
atzer zum
Zeitpunkt
t
2f
1
;:::;T
g
des Regressionsko eÆzienten, b ei dem alle Beobach-
tungen bis auf die Beobachtung (
y
it
;x
it
)indas Mo dell eingehen.
Die Anzahl aller Beobachtungen ist
N
=
N
X
i
=1
n
i
:
Weiter ist
l
(
) eine Verlustfunktion, b eispielsweise werden die folgenden Ver-
lustfunktionen verwendet:
Quadratische Verlustfunktion
l
(
y
it
;
^
(
i;t
)
it
)=(
y
it
^
(
i;t
)
it
)
2
Kullback-Leibler Verlustfunktion
l
(
y
it
;
^
(
i;t
)
it
)=
y
it
ln
(
y
it
^
(
i;t
)
it
)
Minimale Pearson-Residuen
l
(
y
it
;
^
(
i;t
)
it
)=
(
y
it
^
(
i;t
)
it
)
2
v
(^
(
i;t
)
it
)
;
wob ei
v
(^
(
i;t
)
it
) die Varianzfunktion ist. Im bin
aren Mo dell gilt b ei-
spielsweise
v
(^
(
i;t
)
it
)= ^
(
i;t
)
it
(1
^
(
i;t
)
it
) und im Poissonmo dell
v
(^
(
i;t
)
it
)=
^
(
i;t
)
it
.
Da in der verwendeten Datenstruktur (
y
it
;x
it
) korrelierte Daten vorliegen
und deshalb keine Unabh
angigkeit der Residuen angenommen werden kann,
sollte dies bei Kreuzvalidierung ber
ucksichtigt werden. Verfahren f
ur die
Bestimmung der Kreuzvalidierung b ei abh
angigen Beobachtungen werden in
van der Linde (1994) vorgestellt.
39
Kapitel 3
Zuf
allige Eekte und
variierende Ko eÆzienten im
Generalisierten Linearen
Mo dell
Im vorherigen Kapitel wurde die Mo dellierung mit variierenden Ko eÆzien-
ten bzw. mit zuf
alligen Eekten vorgestellt. Beide Ans
atze f
ur Mo dell-
erweiterungen im Generalisierten Linearen Mo dell werden in diesem Kapitel
miteinander kombiniert.
Die Ausgangssituation bei der Verwendung von zeitvariierenden Ko eÆzi-
entenmo dellen ist, da im Zeitverlauf meist die gleichen Beobachtungsein-
heiten o der -grupp en betrachtet werden. Auf Grund der Mewiederholun-
gen o der der vorhandenen Grupp enstruktur, welche b eispielsweise die fa-
mili
are, regionale o der administrative Zugeh
origkeit der einzelnen Beobach-
tungseinheit wiedergibt, ist die Voraussetzung der homogenen Varianzstruk-
tur nicht erf
ullt. Werden diesen Mo dellen zuf
allige Eekte hinzugef
ugt,
40
so wird nicht nur die Grupp enstruktur bzw. die Struktur der Mewieder-
holungen ber
ucksichtigt, sondern dieser Mo dellansatz erm
oglicht auch den
grupp en- bzw. individuensp ezischen Einu der einzelnen Beobachtungs-
einheit zu quantizieren.
Werden die Daten f
ur ein Generalisiertes Lineares Mo dell mit zuf
alligen Ef-
fekten in einem Zeitverlauf erhob en, so kann nichtvorausgesetzt werden, da
der Einu der Ko eÆzienten zeitunabh
angig ist. Um den zeitsp ezischen
Einu der Kovariablen zu ber
ucksichtigen, k
onnen zeitvariierende Ko eÆzi-
enten in das Mo dell aufgenommen werden.
Es ist somit nur nat
urlich, sowohl den Mo dellansatz mit zuf
alligen Eekten
durch die Aufnahme von zeitvariierenden Ko eÆzienten zu erweitern, als auch
b eim Mo dellansatz mit zeitvariierenden Ko eÆzienten zus
atzlichzuf
allige Ef-
fekte zu b etrachten.
41
3.1 Mo dellans
atze und Parametersch
atzungen
Bei der Betrachtung eines Generalisierten Linearen Mo dells, das sowohl zeit-
variierende Ko eÆzienten als auchzuf
allige Eekte b er
ucksichtigt, setzten wir
folgende Datenstruktur voraus.
Die Daten liegen in der Form (
y
it
;x
it
), f
ur
i
= 1
;:::;N
Beobachtungsein-
heiten bzw. -grupp en zu den entsprechenden Zeitpunkten
t
2 f
1
;:::;T
g
vor.
Bei der Mo dellierung wird im folgenden vorausgesetzt, da der Regressions-
ko eÆzient
= (
(1)
;:::;
(
T
)) den zeitvariierenden und die zuf
alligen Ef-
fekte
b
1
;:::;b
N
den individuen-, grupp en- o der auchumweltsp ezischen Ein-
u auf die Zielgr
oe b eschreib en.
3.1.1 Volles Mo dell, parametrische und semi-
parametrische Ans
atze
Ist die Anzahl der b etrachteten Zeitpunkte gering, so kann, wenn keine
zuf
alligen Eekte b er
ucksichtigt werden, separat f
ur jeden b etrachteten Zeit-
punkt ein Mo dell gesch
atzt werden. Dab ei werden immer die Beobachtun-
gen betrachtet, die zu dem zugeh
origen Zeitpunkt erhob en wurden. Unter
Ber
ucksichtigung von zuf
alligen Eekten kann mit diesem Ansatz das folgen-
de Mo dell b etrachtet werden.
Ist
O
die Matrix, deren Ko eÆzienten alle 0 sind, so lautet der lineare Pr
adiktor
im vollen Mo dell
it
=
~
Z
it
+
W
it
b
i
mit
~
Z
it
=(
O
;:::;Z
it
;:::;
O
)und
=(
(1)
;:::;
(
T
)) sowie
b
i
iid
N
(0
;Q
).
42
Bei diesem Mo dellansatz kann die Sch
atzung der Parameter durch den vorge-
stellen EM-Algorithmus mit direkter Gau-Hermite Integration ohne weitere
Anpassungen durchgef
uhrt werden.
Die im vorherigen Kapitel vorgestellen parametrischen und semi-parametrischen
Mo delle lassen sich eb enso einfacherweitern. Der lineare Pr
adiktor
it
=
Z
it
wird durch die Aufnahme der zuf
alligen Eekte durch
it
=
Z
it
(
t
)+
W
it
b
i
mit
b
i
iid
N
(0
;Q
)
ersetzt. Wob ei f
ur
(
t
)eine der in 2.2.1 vorgestellten Parametrisierung
gew
ahlt werden kann.
Mit zunehmender Anzahl der zu b etrachteten Zeitpunkte wird die Anzahl der
zu sch
atzenden Parameter meist so gro, da numerische Sch
atzprobleme b ei
der Bestimmung der Parametersch
atzungen auftreten. Deshalb k
onnen die
ob en vorgestellen Mo delle nur sehr eingeschr
ankt verwendet werden. Der
folgende Lokale Likeliho o dansatz unter Ber
ucksichtigung von zuf
alligen Ef-
fekten ist exibler, da durchdieVerwendung der Kernfunktion mit der Wahl
des Gl
attungsparameters eine b eliebige zeitliche Gl
attung erfolgen kann.
3.1.2 Lokaler Likeliho o dansatz mit zuf
alligen Eekten
Im folgenden stellen wir zwei Algorithmen vor, die sowohl zuf
allige Eekte als
auch zeitvariierende Ko eÆzienten, die mit einem Lokalem Likeliho o dansatz
gesch
atzt werden, ber
ucksichtigen. F
ur die Ber
ucksichtigung der zuf
alligen
Eekte wird zuerst der EM-Algorithmus mit Gau-Hermite Integrationstech-
nik verwendet und im zweiten vorgestellten Algorithmus der EM-Typ e ver-
wendet.
Der lineare Pr
adiktor des Mo dellansatzes lautet in beiden Varianten:
43
Z
it
+
W
it
b
i
mit
b
i
iid
N
(0
;Q
)und
=(
(1)
;:::;
(
T
))
F
ur die zuf
alligen Eekte
b
1
;:::;b
N
,dieden individuen- o der grupp ensp ezi-
schen Einu b eschreib en, wird keine zeitliche Abh
angigkeit vorausgesetzt.
Das heit einerseits, da die zuf
alligen Eekte
b
1
;:::;b
N
unabh
angig vom b e-
trachteten Zeitpunkt
t
2
1
;:::T
sind. Anderseits mu b ei der Sch
atzung des
zeitvariierenden RegessionskoeÆzienten
= (
(1)
;:::;
(
T
))) der Einu
der zuf
alligen Eekte
b
1
;:::;b
N
f
ur jeden Zeitpunkt
t
2
1
;:::;T
gleichblei-
bend ber
ucksichtigt werden.
EM-Algorithmus mit Gau-Hermite Integrationstechnik
In dem zuerst vorgeschlagenen Algorithmus wird ein iteratives Verfahren ver-
wendet, das den zeitvariierenden Ko eÆzienten
sch
atzt, indem die zuf
alligen
Eekte konstant gehalten werden, d. h.
W
it
b
i
gehtals OFFSET in die Mo-
dellgleichung ein. Die Sch
atzung des Parametervektors
=(
(1)
;:::;
(
T
)))
kann mit dem Lokalen Likeliho o dansatz erfolgen.
Anschlieend wird die Kovarianzmatrix
Q
gesch
atzt, wob ei
Z
it
konstant
gehalten wird. Der im vorherigen Kapitel vorgestellte EM-Algorithmus mit
Gau-Hermite Integration kann ohne weitere Anpassungen
ub ernommen wer-
den.
Da nun sowohl
als
Q
bekanntsind,k
onnen die zuf
alligen Eekte
b
1
;:::;b
N
b estimmtwerden.
Dieses Verfahren mu iterativ angewendet werden und es ergibt sich folgender
Ablauf:
Um Startwerte f
ur die Iteration zu erhalten, wird das Mo dell mit zu-
f
alligen Eekten ohne variierenden Ko eÆzienten verwendet.
44
Der lineare Pr
adiktor diese Mo dells hat die Form
it
=
Z
it
+
W
it
b
i
mit
b
i
iid
N
(0
;Q
)
:
{
Zuerst werden
^
und
^
Q
gesch
atzt.
{
Anschlieend erfolgt die Bayes-Sch
atzung der zuf
alligen Eekte
^
b
1
;:::;
^
b
N
.
Unter Verwendung der Startwerte
^
b
h
0
i
i
:=
^
b
i
f
ur
i
= 1
;:::;N
wird die
folgende Iteration f
ur
p
=0
;
1
;:::
durchgef
uhrt:
{
Die zeitvariierenden Ko eÆzienten
^
(
t
)
h
p
+1
i
werden durch den Lo-
kalen Likeliho o dansatz gesch
atzt, bei dem die zuf
alligen Eekte
W
it
^
b
h
p
i
i
konstant gehalten werden. Das heit es gilt
it
=
Z
it
(
t
)
h
p
+1
i
+
OF F S E T
W
it
^
b
h
p
i
i
:
{
Die Aktualisierung der zuf
alligen Eekte
^
b
h
p
+1
i
i
f
ur
i
= 1
;:::;N
erfolgt in zwei Schritten:
Zuerst wird der Sch
atzer
^
Q
h
p
+1
i
durch ein Mo dell mit zuf
alligen
Eekten, b ei dem die zeitvariierenden Ko eÆzienten
^
(
t
)
h
p
+1
i
Teil
des OFFSETS sind, b estimmt. Es wird das Mo dell
it
=
OF F S E T
Z
it
^
(
t
)
h
p
+1
i
+
W
it
^
b
h
p
i
i
mit
b
i
N
(0
;Q
h
p
+1
i
)b etrachtet und anschlieend erfolgt die Ak-
tualisierung der zuf
alligen Eekte
^
b
h
p
+1
i
i
f
ur
i
=1
;:::;N
.
Die Iteration wird b eendet, falls
jj
^
(
t
)
h
p
+1
i
^
(
t
)
h
p
i
jj
jj
^
(
t
)
h
p
i
jj
<
gilt.
45
EM-Typ e Algorithmus
Da beim EM-Typ e Algorithmus eine simultane Sch
atzung des Parameter-
vektors
= (
(1)
;:::;
(
T
)) und der zuf
alligen Eekte
b
i
,
i
= 1
;::: N
mit
b
i
iid
N
(0
;Q
) durchgef
uhrt wird, k
onnen zeitvariierenden Ko eÆzienten pro-
blemlos b er
ucksichtigt werden.
Ist
Q
h
p
i
aus den vorherigen Iterationsschritt b ekannt, so gilt f
ur die Likeliho o d
mit
Æ
0
=(
0
(1)
;:::;
0
(
T
)
;b
1
;:::;b
N
):
l
(
Æ
) =
N
X
i
=1
ln
f
(
y
i
j
b
i
;
)
1
2
N
X
i
=1
b
0
i
Q
h
p
i
1
b
i
=
N
X
i
=1
T
X
s
=1
w
(
t; s
) ln
f
(
y
is
j
b
i
;
(
t
))
1
2
N
X
i
=1
b
0
i
Q
h
p
i
1
b
i
und f
ur die Scorefunktion folgt:
@l
(
Æ
)
@
(
t
)
=
N
X
i
=1
T
X
s
=1
w
(
t; s
)
Z
0
is
D
is
(
(
t
)
;b
i
)
1
is
(
(
t
)
;b
i
) (
y
is
is
(
(
t
)
;b
i
))
@l
(
Æ
)
@b
i
=
T
X
t
=1
W
0
it
D
it
(
(
t
)
;b
i
)
1
it
(
(
t
)
;b
i
) (
y
it
it
(
(
t
)
;b
i
))
Q
h
p
i
1
b
i
Damit hat die Fishermatrix die folgende partitionierte Darstellung:
46
F
(
Æ
)
h
p
i
=
2
6
6
6
6
6
6
6
6
6
6
6
6
6
4
F
1
1
O
F
1
b
1
F
1
b
N
.
.
.
.
.
.
.
.
.
O
F
T
T
F
T
b
1
F
T
b
N
F
b
1
1
F
b
1
T
F
b
1
b
1
.
.
.
.
.
.
.
.
.
O
.
.
.
.
.
.
O
.
.
.
F
b
N
1
F
b
N
T
F
b
n
b
n
3
7
7
7
7
7
7
7
7
7
7
7
7
7
5
mit:
F
t
t
=
N
X
i
=1
T
X
s
=1
w
(
t; s
)
Z
0
is
D
is
(
(
t
)
;b
i
)
1
is
(
(
t
)
;b
i
)
D
0
is
(
(
t
)
;b
i
)
Z
is
F
t
b
i
=
T
X
s
=1
w
(
t; s
)
Z
0
is
D
is
(
(
t
)
;b
i
)
1
is
(
(
t
)
;b
i
)
D
0
is
(
(
t
)
;b
i
)
W
is
F
b
i
t
=
F
0
t
b
i
F
b
i
b
i
=
T
X
t
=1
W
0
it
D
it
(
(
t
)
;b
i
)
1
it
(
(
t
)
;b
i
)
D
0
it
(
(
t
)
;b
i
) +
Q
h
p
i
1
Diese Blo ckgestalt der Fishermatrix kann, wie im Anhang in Abschnitt 6.1.2
gezeigt wird, f
ur eine Vereinfachung bei der Berechnung der Sch
atzer aus-
gen
utzt werden.
47
3.2 Varianz der Parametersch
atzungen
Um inhaltlich die gesch
atzten Parameterwerte interpretieren zu k
onnen, mu
neb en der Sch
atzung der Parameterwerte vor allem die Varianz dieser Sch
atz-
werte ber
ucksichtigt werden. Im folgenden betrachten wir die Sch
atzungen
der Varianzen der Parametersch
atzungen f
ur die ob en vorgeschlagenen Mo-
delle.
Im Generalisierten Linearen Mo dell ohne zuf
allige Eekte und ohne variie-
rende Ko eÆzienten gilt f
ur den gesch
atzten Parametervektor
^
:
^
N
(
; F
1
(
^
))
:
Das heit, die Sch
atzung der Kovarianz kann durch
c
cov
(
^
) =
F
1
(
^
) erfol-
gen, wob ei
F
1
(
^
) die inverse Fishermatrix bezeichnet. Bei diesem Mo del-
lansatz liegen Beobachtungen (
y
i
;x
i
)mit
i
= 1
;:::N
vor und der lineare
Pr
adiktor lautet
i
=
Z
i
,so da f
ur die Fishermatrix
F
(
)=
N
X
i
=1
F
i
(
)=
N
X
i
=1
Z
0
i
D
i
1
i
D
0
i
Z
i
;
mit
D
i
=
D
i
(
)=
@h
(
i
)
@
und
i
=
i
(
)=
cov
(
y
i
j
)gilt.
3.2.1 Varianzsch
atzungen bei Mo dellen mit zuf
alligen
Eekten
Im Generalisierten Linearen Mo dell mit zuf
alligen Eekten sind im linea-
ren Pr
adiktor neb en den gesch
atzten festen Parametern
und
Q
auch die
zuf
alligen Eekte
b
i
f
ur
i
=1
;:::;N
enthalten.
48
EM-Algorithmus mit Gau-Hermite Integrationstechnik
Beim EM-Algorithmus mit Gau-Hermite Integrationstechnik kann wegen
der verwendeten Mischungsdichte f
ur die Varianzsch
atzung nicht die Fisher-
matrix aus dem letzten Iterationsschritt verwendet werden. In Fahrmeir &
Tutz (1994) und in Gourieroux & Montfort (1989) wird vorgeschlagen die
Fishermatrix mit Hilfe der Scorefunktion zu approximieren.
Im Anhang in Abschnitt 6.1.2 wird auf die Parametersch
atzung f
ur den Mo-
dellansatz mit dem linearen Pr
adiktor
it
=
Z
it
+
W
it
b
i
mit
b
i
iid
N
(0
;Q
)
eingegangen. Dieser Mo dellansatz wird mit
b
i
=
Q
1
=
2
a
i
und
=
vec
(
Q
1
=
2
)
sowie
0
=(
0
;
0
) mit
=(
(1)
;:::;
(
T
)) unter Verwendung des Kronecker
Pro duktes
und der Cholesky-Zerlegung von
Q
mit
Q
=
Q
1
=
2
Q
T=
2
umpa-
rametrisiert zu
it
= [
Z
it
; a
0
i
W
it
]
mit
a
i
iid
N
(0
;
I)
:
Ist
M
die Anzahl der St
utzstellen und
c
GH
ij
=
c
GH
ij
(
) die Gewichte der
i
-ten
Beobachtungseinheit f
ur St
utzstelle
j
der Gau-Hermite Integration, so die
Scorefunktion der
i
-ten Beobachtungseinheit mit
s
i
(
)
M
X
j
=1
c
GH
ij
@
ln
f
(
y
i
j
d
j
;
)
@
=
T
X
t
=1
M
X
j
=1
c
GH
ij
@
ln
f
(
y
it
j
d
j
;
)
@
approximiert werden. Die Kovarianz des Parametervektors ^
0
=(
^
0
;
^
0
)kann
somit durch
c
cov
(^
)
c
cov
SC
(^
):=
F
SC
(^
)
1
mit
F
SC
(^
)=
N
X
i
=1
s
i
(^
)
s
0
i
(^
)
49
gesch
atzt werden.
Mit diesem Ansatz l
at sich die Scorefunktion unter Verwendung des linearen
Pr
adiktors
itj
=
Z
it
; d
0
j
W
it
durch
s
GH
i
(
) =
T
X
t
=1
M
X
j
=1
c
GH
ij
(
)
Z
it
; d
0
j
W
it
0
D
itj
1
itj
(
y
it
itj
)
mit
D
itj
=
@h
(
itj
)
@
,
itj
=
h
(
itj
)und
itj
=
cov
(
y
it
j
d
j
;
)b erechnen.
Volles Mo dell
F
ur das volle Mo dell kann der Ansatz ohne Ver
anderungen
ub ernommen werden. Es wurde f
ur die die Sch
atzung von
0
= (
0
;
0
)
mit
= (
(1)
;:::;
(
T
)) ein Generalisiertes Lineares Mo dell mit zuf
alligen
Eekten und dem folgenden linearem Pr
adiktor
it
=
h
(
~
Z
it
+
b
i
)mit
~
Z
it
=(
O
;:::;Z
it
;:::;
O
)
;
verwendet. Die Scorefunktion ist somit approximiert durch
s
i
(
)
s
GH
i
(
)=
M
X
j
=1
c
GH
ij
(
)
@
ln
f
(
y
i
j
d
j
;
)
@
=
T
X
t
=1
M
X
j
=1
c
GH
ij
(
)
h
~
Z
it
; d
0
j
W
it
i
0
D
itj
1
itj
(
y
it
itj
)
und dab ei ist
D
itj
=
@h
(
itj
)
@
itj
=
cov
(
y
it
j
d
j
;
)
itj
=
h
(
itj
)mit
itj
=
h
~
Z
it
; d
0
j
W
it
i
0
B
B
B
B
B
@
(1)
.
.
.
(
T
)
1
C
C
C
C
C
A
:
50
Delta-Metho de
Die Delta-Metho de wird verwendet, um die Sch
atzung
der Varianzen der Parametersch
atzungen im Mo dellen mit parametrischen
und semi-parametrischen zeitvariierenden Ko eÆzienten zu b erechnen.
Wurde b eispielsweise
j
(
t
) parametrisiert als ein Polynom von
t
mit Grad
s
(vgl. Abschnitt 2.2.1), so gilt
j
(
t
) =
P
s
r
=0
jr
t
r
und f
ur die Kovarianz von
cov
(
j
(
t
)) folgt
cov
(
j
(
t
))) =
tcov
(
j
)
t
0
:
Zudem wird die Delta-Metho de auchverwendet, um die Kondenzb
ander f
ur
die Sch
atzung der Varianz
d
var
(^
)mit ^
=
h
(^
)=
h
(
Z
^
) zu b estimmen. Mit
D
(
^
)=
@h
(^
)
@
^
folgt
d
var
(^
)=
D
(
^
)
0
Z
c
cov
(
^
)
Z
0
D
(
^
)
:
Lokaler Likeliho o dansatz
Bei der Parametersch
atzung durch den Lo-
kalen Likeliho o dansatz bei dem die zuf
alligen Eekte durch Gau-Hermite
Integrationstechnik gesch
atzt werden, erfolgt die Sch
atzung von
(
t
)lokal
f
ur jeden b etrachteten Zeitpunkt
t
2
1
;:::T
. Dagegen wird die Sch
atzung
von
global, d. h. unabh
angig von
t
durchgef
uhrt.
Beim Verfahren f
ur die Berechnungen der Varianzen der Parametersch
atzungen
mu deshalb folgende Unterscheidung getroen werden:
Der Sch
atzung von
(
t
)liegt ein Generalisiertes Lineares Mo dell mit
zuf
alligen Eekten zugrunde. Das heit, die Approximation der lokalen
Scorefunktion erfolgt mit
s
i
(
(
t
)
;
)
T
X
s
=1
M
X
j
=1
c
GH
ij
(
(
t
)
;
)
w
(
s; t
)
@
ln
f
(
y
is
j
d
j
;
(
t
)
;
)
@
(
t
)
=
T
X
s
=1
M
X
j
=1
c
GH
ij
(
(
t
)
;
)
w
(
s; t
)
Z
is
; d
0
j
W
is
0
D
isj
1
isj
(
y
is
isj
)
;
51
wob ei
w
(
s; t
) die Gewichte der gew
ahlten Kernfunktion b ezeichnen
und
D
isj
=
@h
(
isj
)
@
isj
=
cov
(
y
is
j
d
j
;
(
t
)
;
)
isj
=
h
(
isj
)mit
isj
=
Z
is
; d
0
j
W
is
(
t
)
!
gilt.
Damit wird die Kovarianzmatrix gesch
atzt mit:
c
cov
(
(
t
)
0
;
0
) =
N
X
i
=1
s
i
(
(
t
)
;
)
s
0
i
(
(
t
)
;
)
!
1
:
Die Sch
atzung von
erfolgt global, so da der folgende marginale An-
satz gew
ahlt werden kann. Mit
= (
(1)
;:::;
(
T
)
;
)und
~
Z
it
=
(
O
;:::;Z
it
;:::;
O
) gilt f
ur die Scorefunktion
s
i
(
)
T
X
t
=1
M
X
j
=1
c
GH
ij
(
)
h
~
Z
it
; d
0
j
W
it
i
0
D
itj
1
itj
(
y
it
itj
)
mit
D
itj
=
@h
(
itj
)
@
itj
=
cov
(
y
it
j
d
j
;
(
t
)
;
)
itj
=
h
(
itj
)mit
itj
=
h
~
Z
it
; d
0
j
W
it
i
0
B
B
B
B
B
@
(1)
.
.
.
(
T
)
1
C
C
C
C
C
A
.
52
und somit folgt f
ur die Kovarianzmatrix
cov
(
) =
N
X
i
=1
s
i
(
)
s
0
i
(
)
!
1
:
Die Kovarianzen
cov
(
(
t
)) ergeb en sich dann aus den entsprechenden Teil-
matrizen von
cov
(
(
t
)
0
;
0
)und die Sch
atzung der Kovarianz von
aus der
entsprechenden Teilmatrix von
cov
(
)=
cov
(
(1)
0
;
;
(
T
)
0
;
0
).
53
Kapitel 4
Simulationsb eispiele
Im vorherigen Kapitel 3 wurden verschiedene Mo dellerweiterungen des Gene-
ralisierten Linearen Mo dells vorgestellt, die gleichzeitig eine zeitliche Variati-
on der Parameter erm
oglichen und grupp en- bzw. individuensp ezische Eek-
te durch die Aufnahme von zuf
alligen Eekten im Mo dellansatz b er
ucksichtigen.
Es wurden verschiedene Sch
atzverfahren f
ur die Bestimmung der Parame-
tersch
atzungen und deren Varianzen angegeb en. Die G
ute der vorgestellten
Verfahren wird in diesem Kapitel durch Simulationsb eispiele veriziert.
Zun
achst werden die vorgestellten Mo dellans
atze f
ur die Sch
atzung der fe-
sten Parameter angewendet und anschlieend wird auf die Bayes-Sch
atzung
der zuf
alligen Eekte eingegangen.
Bei der Sch
atzung der festen Parameter werden der Lokale Likeliho o dansatz
und der Ansatz des vollen Mo dells miteinander verglichen. Dab ei zeigt sich,
da der Lokale Likeliho o dansatz, b ei dem f
ur die Sch
atzung der zuf
alligen Ef-
fekte die Gau-Hermite Integrationstechnik verwendet wird, sehr gute Resul-
tate liefert. Wird dagegen der EM-Typ e Algorithmus angewendet, so liefert
dieser nur zufriedenstellende Ergebnisse, wenn die Streuung der zuf
alligen Ef-
54
fekte sehr gering ist. Damit ist eine Aussage
ub er die Sch
atzung der zuf
alligen
Eekte b eim EM-Typ e Algorithmus nur b edingt aussagekr
aftig.
55
4.1 Sch
atzung der festen Parameter
In den folgenden Simulationstudien werden die Ergebnisse aus dem vollen
Mo dell mit denen des Lokalen Likeliho o dansatzes mit Gau-Hermit Integra-
tionstechnik bzw. mit EM-Typ e Algorithmus jeweils unter Ber
ucksichtigung
von zuf
alligen Eekten vorgestellt.
Aufbau der Simulationsstudie
In den Simulationsstudien wird das bin
are logistische Mo dell
it
=
h
(
it
)=
exp(
)
1 + exp (
)
mit
it
=
0
(
t
)+
1
(
t
)
x
it
+
b
i
und
b
i
N
(0
;Q
) mit 100 Beobachtungseinhei-
ten f
ur 10 Zeitpunkte gew
ahlt.
Der Regressor
x
it
ist mit
x
it
=
x
i
f
ur
t
= 1
;:::;
10 und
x
i
N
(0
;
1) f
ur
i
=1
;:::;
100 gew
ahlt worden. Die Parametervorgab en f
ur
sind
0
= (
0
(1)
;
0
(2)
;:::;
0
(10))
0
=(1
;
0
:
9
;:::;
0
:
1)
0
1
= (
1
(1)
;
1
(2)
;:::;
1
(10))
0
=(0
:
1
;
0
:
2
;:::;
1)
0
:
F
ur
wird b eim vollen Mo dell und b eim Lokalen Likeliho o dansatz mit Gau-
Hermite Integrationstechnik die Parametervorgab e mit
=
p
Q
= 2, dagegen
b eim Lokalen Likeliho o dansatz unter Verwendung des EM-Typ e Algorithmus
mit
=
p
Q
=
p
(0
:
1) gew
ahlt.
4.1.1 Das volle Mo dell mit zuf
alligen Eekten
Das volle Mo dell unter Ber
ucksichtigung von zuf
alligen Eekten wurde im
Abschnitt 3.1.1 vorgestellt. Die Sch
atzung der Parameter
und
=(
0
(1)
;
56
1
(1)
;:::;
0
(10)
;
1
(10)) erfolgt durch den Mo dellansatz eines Generalisier-
ten Linearen Mo dells mit zuf
alligen Eekten und dem linearem Pr
adiktor
it
=
h
(
~
Z
it
+
b
i
)mit
~
Z
it
=(
O
;:::;Z
it
;:::;
O
)mit
Z
it
=(1
;x
it
)
:
Sch
atzung der Parametervarianzen
Die Sch
atzung der Varianz wird, wie in Abschnitt 3.2.1 vorgestellt, durchdie
Approximation der inversen Fishermatrix mit der Scorefunktion, d. h.
c
cov
(^
)
F
SC
1
mit
F
SC
=
n
X
i
=1
s
i
(^
)
s
0
i
(^
)
durchgef
uhrt. Ist
M
die Anzahl der St
utzstellen f
ur die numerische Integra-
tion so gilt:
s
i
(
)
s
GH
i
(
)=
M
X
j
=1
c
GH
ij
(
)
@
ln
f
(
y
i
j
d
j
;
)
@
=
T
X
t
=1
M
X
j
=1
c
GH
it
(
)
~
Z
it
; d
j
0
D
itj
(
itj
)
1
it
(
itj
) (
y
it
it
(
itj
))
mit
itj
=
~
Z
it
; d
j
In den folgenden Abbildungen werden die Ergebnisse der Simulationsstudie
graphisch dargestellt. Dab ei wird deutlich, da die empirischen Standardab-
weichungen
St
emp
(
(
t
)) :=
v
u
u
t
1
U
1
U
X
u
=1
^
u
(
t
)
(
t
)
2
bzw.
St
emp
(
) :=
v
u
u
t
1
U
1
U
X
u
=1
^
u
2
57
sehr gut mit den entsprechenden gesch
atzten Standardabweichungen
St
SC
(
(
t
)) :=
1
U
U
X
u
=1
r
d
var
F
SC
^
u
(
t
)
bzw.
St
SC
(
) :=
1
U
U
X
u
=1
r
d
var
F
SC
^
u
ub ereinstimmen. Mit
(
t
) =
1
U
U
P
u
=1
^
u
(
t
) und
(
t
) =
1
U
U
P
u
=1
^
u
(
t
) werden die
durchschnittlichen Parametersch
atzungen b ezeichnet.
Ergebnisse der Simulationsstudie
In den Abbildungen 4.1 und 4.2 werden f
ur die Parameter
0
und
1
zu jedem
betrachteten Zeitpunkt
t
= 1
;
2
;:::;
10 alle Parametersch
atzungen in Form
eines Box-Plots eingetragen. Zus
atzlich werden die Mo dellvorgab en sowie
die empirischen und gesch
atzten Standardabweichungen mit
Modellvorgabe
+− empirische Standardabweichung
+− geschätzte Standardabweichung
gekennzeichnet.
58
246810
-0.5 0.0 0.5 1.0 1.5 2.0
marginal model
beta 0
Anzahl der Simulationen
U= 160
Abbildung 4.1: Simulationsergebnisse f
ur
0
(1)
;
0
(2)
;:::;
0
(10) im vollen
Mo dell
246810
-1 0 1 2
marginal model
beta 1
Anzahl der Simulationen
U= 160
Abbildung 4.2: Simulationsergebnisse f
ur
1
(1)
;
1
(2)
;:::;
1
(10) im vollen
Mo dell
59
F
ur die zugeh
origen Parametersch
atzungen von
wurden die folgenden durch-
schnittlichen Parameter- bzw. Varianzwerte berechnet.
=
1
U
U
P
u
=1
^
u
St
emp
(
)
St
SC
(
)
2 2.024 0.227 0.214
4.1.2 Lokaler Likeliho o dansatz mit zuf
alligen Eekten
und Gau-Hermite Integrationstechnik
In diesem Abschnitt werden die Ergebnisse der Parametersch
atzungen beim
Lokalen Likeliho o dansatz unter Ber
ucksichtigung von zuf
alligen Eekten vor-
gestellt, wob ei die Gau-Hermite Integrationstechnik verwendet wird.
Die Sch
atzung von
0
:= (
0
;
)und
= (
0
(1)
;
1
(1)
;:::;
0
(10)
;
1
(10))
erfolgt somit durch ein Generalisiertes Lineares Mo dell mit zuf
alligen Ef-
fekten und (zeit-)variierenden Ko eÆzienten unter Verwendung des Lokalen
Likeliho o d Ansatzes mit dem linearen Pr
adiktor:
it
=
h
(
Z
it
(
t
)+
b
i
)mit
Z
it
=(1
;x
it
),
(
t
)=(
0
(
t
)
;
1
(
t
)) und
b
i
N
(0
;Q
)
:
Es wird das in Abschnitt 3.1.2 vorgestellte Verfahren verwendet, wob ei f
ur
die Kernfunktion der Gau-Kern mit
=1 und als Startwert
h
0
i
=2 sowie
f
ur das Abbruchkriterium der Iteration
=0
:
05 gew
ahlt wurden.
Sch
atzung der Parametervarianzen
Bei der Parametersch
atzung wird der Ko eÆzienten
(
t
) mit Hilfe der lokalen
Likeliho o d und den Gewichten des Gaukerns berechnet, dagegen erfolgte
die Parametersch
atzung von
global, d. h. ohne lokale Gewichtung. Des-
halb mu b ei der Sch
atzung der Varianzen f
ur die Parametersch
atzungen f
ur
60
var
(
(
t
)) und
var
(
)unterschiedliche Verfahren angewendet werden, die in
Abschnitt 3.2.1 vorgestellt wurden.
Die Lokale Sch
atzung von
var
(
(
t
)) erfolgt mit:
cov
(
(
t
)
;
) =
N
X
i
=1
s
i
(
(
t
)
;
)
s
0
i
(
(
t
)
;
)
!
1
und
s
i
(
(
t
)
;
)
T
X
s
=1
M
X
j
=1
c
GH
ij
(
(
t
)
;
)
w
(
s; t
)
@
ln
f
(
y
is
j
d
j
;
(
t
)
;
)
@
(
t
)
=
T
X
s
=1
M
X
j
=1
c
GH
ij
(
(
t
)
;
)
w
(
s; t
)
Z
is
; d
0
j
W
is
0
D
isj
1
isj
(
y
is
isj
)
:
F
ur die globale Sch
atzung von
mit
=(
(1)
;:::;
(
T
)
;
)und
~
Z
it
=
(
O
;:::;Z
it
;:::;
O
) gilt:
cov
(
) =
N
X
i
=1
s
i
(
)
s
0
i
(
)
!
1
mit
s
i
(
)
T
X
t
=1
M
X
j
=1
c
GH
ij
(
)
h
~
Z
it
; d
0
j
W
it
i
0
D
itj
1
itj
(
y
it
itj
)
:
F
ur die Sch
atzung der Kovarianz von
(
t
) bzw.
m
ussen die entsprechenden
Teilmatrizen von
cov
(
) b etrachtet werden.
Ergebnisse der Simulationsstudie
Die Ergebnisse der Simulation werden f
ur jeden b etrachteten Zeitpunkt
t
=
1
;
2
;:::;
10 in Abbildung 4.3 und 4.4 in Form eines Box-Plots eingetragen.
61
Mit
Modellvorgabe
+− empirische Standardabweichung
+− geschätzte Standardabweichung
werden zus
atzlich die Mo dellvorgab en sowie die empirischen und gesch
atzten
Standardabweichungen gekennzeichnet.
246810
0.0 0.5 1.0
local Likelihood
beta 0
Anzahl der Simulationen
U= 100
Abbildung 4.3: Simulationsergebnisse f
ur
0
(1)
;
0
(2)
;:::;
0
(10) b eim Loka-
len Likeliho o dansatz mit zuf
alligen Eekten und Gau-Hermite Integrations-
technik und
=1
62
246810
−0.5 0.0 0.5 1.0 1.5
local Likelihood
beta 1
Anzahl der Simulationen
U= 100
Abbildung 4.4: Simulationsergebnisse f
ur
1
(1)
;
1
(2)
;:::;
1
(10) b eim Loka-
len Likeliho o dansatz mit zuf
alligen Eekten und Gau-Hermite Integrations-
technik und
=1
Die Ergebnisse der Parametersch
atzungen von
werden in der folgenden
Tab elle zusammengestellt:
=
1
U
U
P
u
=1
^
u
s
1
U
1
U
P
u
=1
^
u
2
1
U
U
P
u
=1
r
d
var
F
SC
^
u
2 1.890 0.199 0.202
Die Struktur der zeitvariierenden Ko eÆzienten wird bei der Sch
atzung der
Parameterwerte exakt wiedergegb en. Da die Mo dellvorgab en im Mittel
nicht genau getroen werden, ist nicht verwunderlich, da beim verwendeten
Algorithmus alle Beobachtungen gewichtet in Abh
angigkeit der verwendeten
Kernfunktion eingehen. Insb esondere sind deshalb zu Beginn (d. h. f
ur
63
t
= 1 und
t
= 2) und zum Ende (d. h. f
ur
t
=9 und
t
= 10) die gesch
atzten
Parametervarianzen gr
oer als in der Mitte des Beobachtungszeitraumes f
ur
t
=4,
t
=5 und
t
=6, bei denen die gesch
atzten und emirischen Varianzen
sehr gut
ub ereinstimmen.
4.1.3 Lokaler Likeliho o dansatz mit zuf
alligen Eekten
und EM-Typ e Algorithmus
Der Lokale Likeliho o dansatz unter Ber
ucksichtigung von zuf
alligen Eekten
mit Anwendung des EM-Typ e Algorithmus wurde in Abschnitt 2.1.2 vorge-
stellt. Bei den durchgef
uhrten Simulationsstudien hat sich gezeigt, da der
EM-Typ e Algorithmus ein sehr schnelles Verfahren ist, und recht gute Er-
gebnisse liefert, wenn die Standardabweichung der normalverteilten zuf
alligen
Eekte sehr klein ist (vgl. Hennevogl, 1991).
Es werden jedo ch keine guten Sch
atzresultate erzielt, wenn die Standard-
abweichung der zuf
alligen Eekte nicht sehr klein ist. Deshalb wurden die
Parametervorgab en f
ur
0
und
1
mit
0
= (
0
(1)
;
0
(2)
;:::;
0
(10))
0
=(1
;
0
:
9
;:::;
0
:
1)
0
1
= (
1
(1)
;
1
(2)
;:::;
1
(10))
0
=(0
:
1
;
0
:
2
;:::;
1)
0
identisch zu den Vorgab en der letzten b eiden Simulationsstudien, jedo chdie
Standardabweichung der normalverteilten zuf
alligen Eekte mit
=
p
Q
=
p
0
:
1, statt mit
= 2 gew
ahlt.
Die Sch
atzung von
0
:= (
0
;
)und
= (
0
(1)
;
1
(1)
;:::;
0
(10)
;
1
(10))
erfolgt durch ein Generalisiertes Lineares Mo dell mit zuf
alligen Eekten und
(zeit-)variierenden Ko eÆzienten unter Verwendung der lokalen Likeliho o d
und dem linearen Pr
adiktor:
it
=
h
(
Z
it
(
t
)+
b
i
)mit
Z
it
=(1
;x
it
),
(
t
)=(
0
(
t
)
;
1
(
t
)) und
b
i
N
(0
;Q
)
64
Ergebnisse der Simulationsstudie
Die Sch
atzungen der Parameter
0
(1)
;
0
(2)
;:::;
0
(10) und
1
(1)
;
1
(2)
; :::;
1
(10) werden in den Abbildungen 4.5 und 4.6 mit Box-Plots dargestellt. So-
wie in Abbildung 4.7 sind die Sch
atzungen f
ur
dargestellt. Die Sch
atzungen
f
ur
sind nicht zufriedenstellend gesch
atzt worden. Ist die Varianz der
zuf
alligen Eekte gr
oer, so werden die Sch
atzungen f
ur den Mo dellpara-
meter
sehr ungenau. Deshalb werden wir dieses Verfahren nicht bei den
Anwendungsb eispielen im anschlieenden Kapitel 5 verwenden.
246810
-0.2 0.0 0.2 0.4 0.6 0.8 1.0 1.2
local Likelihood mit EM-Type
beta 0
gamma = 1
Anzahl der Simulationen
U= 50
Abbildung 4.5: Simulationsergebnisse f
ur
0
(1)
;
0
(2)
;:::;
0
(10) b eim Loka-
len Likeliho o dansatz mit zuf
alligen Eekten und EM-Typ e Algorithmus
65
246810
0.0 0.2 0.4 0.6 0.8 1.0 1.2
local Likelihood mit EM-Type
beta 1
gamma = 1
Anzahl der Simulationen
U= 50
Abbildung 4.6: Simulationsergebnisse f
ur
1
(1)
;
1
(2)
;:::;
1
(10) b eim Loka-
len Likeliho o dansatz mit zuf
alligen Eekten und EM-Typ e Algorithmus
0.0 0.5 1.0 1.5 2.0
0.30 0.35 0.40 0.45 0.50
Anzahl der Simulationen
U= 50
local Likelihood mit EM-Type
theta
gamma = 1
Startwert = 0.5
Abbildung 4.7: Simulationsergebnisse f
ur
beim Lokalen Likeliho o dansatz
mit zuf
alligen Eekten und EM-Typ e Algorithmus
66
4.2 Bayes-Sch
atzung der zuf
alligen Eekte
Die Bayes-Sch
atzung der zuf
alligen Eekte kann mit dem in Abschnitt 2.1.1
vorgestellten Verfahren durchgef
uhrt werden.
In Abbildung 4.8 sind die Simulationsergebnisse f
ur die Bayes-Sch
atzung der
zuf
alligen Eekte
b
i
dargestellt. Dab ei ist ein bin
ares logistisches Mo dell mit
zuf
alligen Eekten und dem Mo dellansatz
E
(
y
it
j
x
it
)=
exp
(
it
)
1+
exp
(
it
)
und
it
=
0
+
x
(1)
it
1
+
x
(2)
it
2
+
b
i
:
gew
ahlt worden.
Bei der Simulation ist die Anzahl der Untersuchungseinheiten
N
= 100. Die
Anzahl der Mewiederholungen ist f
ur alle Untersuchungseinheiten konstant
mit
n
i
=4 festgelegt worden.
Der Regressor
x
(1)
i
ist bin
ar mit
x
(1)
i
2 f
0
;
1
g
und
x
(2)
i
ist metrisch durch
x
(2)
i
N
(0
;
10) gew
ahlt worden. F
ur alle Mewiederholungen
t
= 1
;:::;
4
sind diese Regressoren nicht ver
andert worden. Damit gilt
(1
;x
(1)
i
;x
(2)
i
):=(1
;x
(1)
it
;x
(2)
it
)
:
Weiter gilt f
ur die zuf
alligen Eekte
b
i
N
(0
;Q
)mit
Q
=
2
=4.
Die Sch
atzungen dieser Eekte f
ur 50 Simulationsdurchl
aufe sind in der
obigen Abbildung 4.8 dargestellt. In vier F
allen konnte selbst nach 100
Iterationsschritten kein
b
i
gefunden werden, welches das Abbruchkriterium
j
b
(
neu
)
i
b
(
alt
)
i
j
j
b
(
neu
)
i
j
0
:
001 erf
ullt. In allen anderen F
allen wurde der Sch
atzwert
nach maximal 8 Iterationen gefunden.
Es ist auallend, da sehr groe bzw. sehr kleine Simulationsvorgab en von
b
i
nicht ausreichend gesch
atzt werden. Der QQ-Plot der Simulationsvorgab en
f
ur
b
i
gegen die Quantile der Normalverteilung
N
(0
;
4) in Abbildung 4.9
67
Personen 1 bis 50
0 1020304050
−4−2024
Simulationsvorgabe
Schätzung
keine Schätzung möglich
Schätzung der Randomeffekte
Personen 51 bis 100
0 1020304050
−4−2024
Abbildung 4.8: Simulationsergebnisse f
ur
zeigt, da die Mo dellvoraussetzung
b
i
N
(0
;
4) in der konkret vorliegenden
Stichprob e nicht exakt gegeb en ist.
Dies ist auch der Grund weshalb die Sch
atzung der festen Parameter nicht
exakt mit den Mo dellvorgab en
ub ereinstimmen. Dies wird in Abbildung
4.10 verdeutlicht, in dem die Kerndichtesch
atzungen f
ur die festen Parameter
0
,
1
und
2
sowie die aufsteigend sortierten Parametersch
atzungen f
ur
dargestellt werden.
68
Quantile der Normalverteilung N(0,4)
Simulationsvorgabe der Randomeffekte
−50 5
−6−4−20 2 4 6
QQ−Plot
Abbildung 4.9: QQ-Plot der Simulationsvorgab en
Kerdichteschätzung
−2.0 −1.0 0.0
0.0 1.5 3.0
beta 0
Kerdichteschätzung
−1012
0.0 1.0 2.0
beta 1
Kerdichteschätzung
−0.6 −0.4 −0.2
0 5 10 20
beta 2
sort. Parameterschätzung
0 1020304050
1.5 1.7 1.9
················································
sigma
Abbildung 4.10: Sch
atzungen der festen Parameter
69
Kapitel 5
Anwendungsb eispiele
Es werden zun
achst zwei aus dem medizinischen Umfeld stammende Anwen-
dungssb eipiele f
ur Poissonmo delle vorgestellt. Anschlieend wird auf zwei
Beispiele f
ur Verweildauermo delle aus dem sozial-
okonomischen bzw. aus
dem medizinischen Bereich eingegangen.
In allen Studien werden die relevanten Merkmale der gleichen Personen zu
verschiedenen Zeitpunkten untersucht, so da es sich aus statistischer Sicht
um Mewiederholungen an den gleichen Beobachtungseinheiten im Zeitver-
lauf handelt. Durch die Verwendung von variierenden Ko eÆzienten kann
dieser zeitliche Einu untersucht werden, und durch die Aufnahme der
zuf
alligen Eekte wird sowohl die Datenstruktur der Mewiederholungen als
auch der individuensp ezische Einu b er
ucksichtigt.
70
5.1 Poissonmo delle
5.1.1 Emesis Datensatz
Datenbeschreibung
H
auge Neb enwirkungen der Chemotherapie sind
Ub elkeit und Erbrechen
(Emesis), die nachAngab en der Patienten sub jektiv als die am meisten be-
lastende Neb enwirkung empfunden wird.
Durch eine neue Art von Medikamenten sollen diese Neb enwirkungen verrin-
gert werden. Die 122 Patientinnen der vorliegenden Studie (Dietz & B
ohning,
1994) mit gyn
akologischen Tumoren wurden alle mit einer Chemotherapie in
Kombination mit einem anti-emetischen Wirksto b ehandelt. Es wurden
dab ei zwei unterschiedliche Wirstoe verwendet. Mit Medikament1(
Med 1
)
wurden 70 und mit Medikament2(
Med 2
) wurden 52 Patientinnen b ehandelt.
Das Empnden der
Ub elkeit einer Patientin ob jektiv zu messen ist sehr
schwierig. Einfacher ist, es die Anzahl der Erbrechungsanf
alle pro Tag zu
z
ahlen, wob ei wiederholtes Erbrechen in einer kurzen Zeitspanne als ein
Anfall gewertet wird. Der Zusammenhang zwischen der Anzahl der Erbre-
chungsanf
alle und dem sub jektiven Empnden der
Ub elkeit ist oensichtlich
gegeb en, so da in dem hier vorgeschlagenen Mo dell die Zielgr
oe die Anzahl
der Erbrechungsanf
alle pro Tag ist.
Es ist sinnvoll, die Patientinnen
ub er mehrere Tage zu b eobachten, um sicher-
zustellen, da ein neues Anti-Emetikum eventuell nichtnur zu einer zeitlichen
Verz
ogerung der Emesis f
uhrt. Die Patientinnen in dieser Studie wurden nach
der Anzahl der Erbrechungsanf
alle pro Tag am Therapietag selbst (
t
= 1)
und an den n
achsten f
unf folgenden Tagen (
t
=2
;:::;t
= 6) b efragt.
Als weitere prognostische Faktoren wurde das Alter der Patientinnen (
age
)
und die Information, ob dies die erste Chemotherapie (
Zyk 1
=1) oder eine
71
weitere Chemotherapie (
Zyk 1
= 0) ist, b er
ucksichtigt.
Mo dellansatz
Die Zielgr
oe, d. h. die Anzahl der t
aglichen Erbrechungsanf
alle, wurde als
eine p oissonverteilte Zufallsvariable gew
ahlt. Es wurde der Mo dellansatz
E
y
it
j
x
Med 1
;x
Med 2
;x
age
;x
Zyk 1
= exp
x
Med 1
it
+
x
Med 2
it
+ln(
x
age
it
)+
x
Zyk 1
it
+
b
i
mit
b
i
iid
N
(0
;
2
)verwendet. Der lineare Pr
adiktor des Mo dells enth
alt kei-
nen Interzept, damit der Einu der unterschiedlichen Medikamente, d. h.
die Parameterwerte von (
Med 1
), bzw. (
Med 2
) direkt interpretiert werden
k
onnen.
Das verwendete Mo dell ist ein Generalisiertes Lineares Mo dell mit zuf
alligen
Eekten und zeitvariierenden Ko eÆzienten. Als Mo dellansatz wird der zen-
trale Ansatz dieser Arb eit, der Lokale Likeliho o dansatz unter Ber
ucksichti-
gung von zuf
alligen Eekten verwendet, wob ei die Berechnung mit Hilfe der
Gau-Hermite Integrationstechnik durchgef
uhrt wird.
In diesem Anwendungsb eispiel werden wir auf die Wahl des Gl
attungspa-
rameters der Kernfunktion eingehen. Dab ei werden verschiedene Optima-
lit
atskriterien b ez
uglich der Wahl des Gl
attungsparameters angewendet und
f
ur unterschiedliche Gl
attungsparameters die Ergebnisse gegen
ub ergestellt.
Ergebnisse
Die Aufgab enstellung der Studie ist die Unterschiede der b eiden Medikamen-
te hinsichtlich der Emetogenit
at zu b eurteilen. Wie die folgende Darstellung
5.1 zeigt, ist die Anzahl der Erbrechungsanf
alle im gesammten Zeitverlauf
b ei der Einnahme von Medikament2 geringer als b ei Medikament 1.
Neb en dem Zeitverlauf der Parametersch
atzungen, werden auch die Kurven
der punktweisen Fehlerb
ander gezeigt. Letzter ergeb en sichaus der Addi-
72
tion, bzw. Subtraktion der Parametersch
atzwerte und den entsprechenden
Varianzen der Parametersch
atzungen.
Times
-2 0 2 4
t = 1 t = 2 t = 3 t = 4 t = 5 t = 6
Med 1
Med 2
gamma = 1
Abbildung 5.1: VergleichderParametersch
atzungen von
Med 1
und
Med 2
Es wurde der Gl
attungsparameter mit
=1, d. h. der optimale Gl
attungs-
parameter nach dem Kreuzvalidierungskriterium mit Kullback-Leibler Ver-
lustfunktion, gew
ahlt.
Bevor wir auf den Einu des Alters und des Therapiezyklusses eingehen,
stellen wir die Ergebnisse der Kreuzvalidierung f
ur die Wahl des Gl
attungs-
parameters vor.
Wahl und Einu des Gl
attungsparameters
Die Wahl des Gl
attungsparameters wurde mit Kreuzvalidierung, die in Ab-
schnitt 2.2.2 bei der Einf
uhrung des Lokalen Likeliho o dansatzes vorgestellt
wurde, durchgef
uhrt. Der optimale Gl
attungsparameter ist der Parameter,
der im Minimum der Verlustfunktion liegt. Wie in Abbildung 5.2 zu sehen
ist, wird der optimale Gl
attungsparameter durch Kreuzvalidierung bei Ver-
wendung der Quadratischen bzw. der Poisson-Verlustfunktion mit
= 1
:
5,
dagegen bei der Kullback-Leibler Verlustfunktion mit
=1 b estimmt.
73
Regressoren:
Med 1, Med 2, Age, Zyk 1
Emesis Datensatz
•
••••••
•
•
•
•
gamma
CV(gamma)
12345
3.0 3.1 3.2 3.3
Quadratische Verlustfunktion
•
•••••••
•
•
•
gamma
CV(gamma)
12345
0.115 0.125 0.135 0.145
Kullback-Leibler Verlustfunktion
•
•
•••••••••
gamma
CV(gamma)
12345
1.2 1.6 2.0 2.4
Pearson Verlustfunktion
Abbildung 5.2: Wahl des Gl
attungsparameters durch Kreuzvalidierung
Durch die Bestimmung des optimalen Gl
attungsparameters mit Kreuzvali-
dierung durch eine der ob en genannten Verlustfunktionen wird sichergestellt,
da die Struktur der zeitabh
angigen Ko eÆzienten erkennbar ist. Es wird
jedo ch
ub er starke zeitliche Variation in den Daten gegl
attet, damit der zeit-
liche Trend optimal zu erkennen ist.
Wird ein Parameterwert f
ur den Gl
attungsparameter gew
ahlt, der deutlich
gr
oer als der optimale Gl
attungsparameter ist, ist die relevante zeitabh
angige
Struktur nicht mehr erkennbar. W
ahlt man dagegen einen Parameterwert
f
ur den Gl
attungsparameter, der kleiner als der optimale Gl
attungsparameter
74
ist, l
at sich die den Daten zugrundeliegende Struktur genauer analysieren.
Aus diesem Grund werden im folgenden die Parametersch
atzungen mit dem
Gl
attungsparameter
= 0
:
5 und
= 1
:
5 in den Abbildung 5.3 und 5.4
dargestellt.
kein optimales gamma
gamma = 0.5
mit Konfidenzband
Emesis Datensatz
Times
-20246
t = 1 t = 2 t = 3 t = 4 t = 5 t = 6
Med
Med 1
Med 2
Times
-1.5 -1.0 -0.5 0.0 0.5
t = 1 t = 2 t = 3 t = 4 t = 5 t = 6
Age
Times
-0.4 -0.2 0.0 0.2
t = 1 t = 2 t = 3 t = 4 t = 5 t = 6
Zyk 1
Abbildung 5.3: Parametersch
atzungen mit
=0
:
5
In dieser Studie hatten Patientinnen mit der ersten Chemotherapieb ehand-
lung weniger Erbrechungsanf
alle als Patientinnen mit zwei oder mehr Be-
handlungen.
Am zweiten Tag nach der Behandlung (
t
= 3) war die Anzahl der Anf
alle
jedo ch gleich. Da die ob ere Grenze des Kondenzbandes auch an den ande-
75
optimales gamma
mit Kreuzvalidierung:
gamma = 1.5
mit Konfidenzband
Emesis Datensatz
Times
-20246
t = 1 t = 2 t = 3 t = 4 t = 5 t = 6
Med
Med 1
Med 2
Times
-1.0 -0.5 0.0 0.5
t = 1 t = 2 t = 3 t = 4 t = 5 t = 6
Age
Times
-0.2 -0.1 0.0 0.1 0.2
t = 1 t = 2 t = 3 t = 4 t = 5 t = 6
Zyk 1
Abbildung 5.4: Parametersch
atzungen mit
=1
:
5
ren Tagen ob erhalb der Null-Linie liegt, ist an allen Tagen der Einu des
Therapiezyklusses nicht signikantvon Null verschieden.
Bei j
ungeren Patientinnen waren die Erbrechungsanf
alle h
auger. Auch die-
ser Eekt ist nicht signikant, da bei dem optimalen Gamma, d. h.
=1
:
5
die Null-Linie
ub erhalb des Kondenzbandes liegt.
Die punktweisen Kondenzintervalle f
ur
=1
:
5 sind zum Teil an den R
andern
untypisch schm
aler, als in der Mitte des Beobachtungszeitraumes.
Betrachtet man die Ergebnisse der Parametersch
atzungen mit
= 0
:
5, so
76
erkennt man, da
ub er eine groe \zackenartige" Variation der Parameter-
sch
atzungen gegl
attet wurde und deshalb die Kondenzintervalle sehr breit
sind.
5.1.2 Epileptic Datensatz
Datenbeschreibung
In Thall & Vail (1990) und Breslow&Clayton (1993) sowie in Diggle, Liang
& Zeger (1994) werden die Ergebnisse aus eine klinischen Studie vorgestellt,
bei der 59 Patienten, die an epileptischen Anf
allen leiden, zus
atzlich zur
Standardb ehandlung entweder ein neues Medikamentoder ein Placeb o ver-
abreichtwurde. Die Patienten hab en sich ab Beginn der Behandlung viermal,
jeweils nach zwei Wochen, in der Klinik vorgestellt.
Die Daten liegen in der Form (
y
it
;x
it
)mit
i
= 1
;:::;
59 (Personen)
t
=
1
;:::;
4 (Mewiederholungen) vor. Dab ei ist
y
it
die Anzahl der epileptischen
Anf
alle der i-ten Person in einem zweiw
ochigen Intervall vor dem t-ten Arzt-
b esuch. F
ur den Regressorenvektor
x
it
=
x
i
= (1
;x
Treat
i
;x
B ase
i
;x
Ag e
i
)der
i-ten Person gilt
x
T reat
=0, wenn keine zus
atzliche Behandlung (Placeb o)
bzw.
x
T reat
= 1, wenn die Behandlung mit dem neuen Medikament erfolgt
ist. Der BASELINE
x
B ase
stellt die patientensp ezische Anzahl der epilepti-
schen Anf
alle dar und ist die logarithmierte Anzahl der epileptischen Anf
alle
vor Beginn der neuen Behandlungsmetho de innerhalb eines Zeitraumes vom
8Wochen.
Mo dellans
atze
Im folgenden werden die Ergebnisse verschiedener Mo dellans
atze vorgestellt,
wob ei jeweils ein Poissonmo dell verwendet wurde. In den Tab ellen werden
77
neb en den Parametersch
atzungen auch die zugeh
origen Standardabweichun-
gen angegeb en.
Beim ersten Mo dellansatz werden keine variierenden Ko eÆzienten und keine
zuf
allige Eekte verwendet.
0
T reat
Base
Ag e
-3.85 (0.39) -0.04 (0.05) 1.22 (0.03) 0.52 (0.1)
Der obige Mo dellansatz kann durch Ber
ucksichtigung von zuf
alligen Eekten
erweitert werden und man erh
alt die folgenden Parametersch
atzungen.
0
Treat
Base
Ag e
-4.54 (0.39) -0.25 (0.04) 1.28 (0.03) 0.74 (0.09) 0.6735 (0.031)
Im n
achsten Schritt werden neb en der Aufnahme von zuf
alligen Eekten
auch zeitvariierende Ko eÆzienten b etrachtet. Dab ei werden drei verschie-
dene Ans
atze f
ur die zeitvariierenden Ko eÆzienten verwendet. Beim ersten
Ansatz wird der Marginale Likeliho o dansatz verwendet, beim zweiten die
Polynomiale Regression mit einem Polynom 3. Grades. Bei diesem Ansatz
k
onnen keine Standardabweichungen f
ur die Parametersch
atzungen angege-
ben werden, da die approximierte Fishermatrix aus numerischen Gr
unden
nicht invertiert werden kann. Beim zuletzt b etrachteten Mo dell wird der
Lokale Likeliho o dansatz mit
=0
:
1verwendet.
Obwohl die drei Mo dellans
atze aufgrund der unterschiedlichen Ber
ucksich-
tigung des zeitlichen Einusses zu verschiedenen Parametersch
atzungen f
ur
die zeitabh
angigen Ko eÆzienten f
uhren, lassen sichdie gleichen Tendenzen
erkennen. Damit ist die qualitative Aussage der Parametersch
atzungen in
78
allen drei Mo dellen gleich. Insb esondere wird der zeitunabh
angige individu-
ensp ezische Einu der zuf
alligen Eekte in allen Mo dellen gleichgesch
atzt.
Volles Mo dell
0
t
T reat
t
B ase
t
Ag e
t
t=1 -6.13 (0.66) -0.23 (0.1) 1.40 (0.05) 1.1 (0.16)
t=2 -2.79 (1.19) -0.19 (0.11) 1.13 (0.09) 0.38 (0.34) 0.67 (0.04)
t=3 -4.99 (0.63) -0.25 (0.09) 1.33 (0.10) 0.83 (0.15)
t=4 -4.03 (0.94) -0.36 (0.16) 1.25 (0.08) 0.61 (0.26)
Polynomiale Regression
0
t
T reat
t
B ase
t
Ag e
t
t=1 -6.13 -0.23 1.40 1.1
t=2 -2.79 -0.19 1.13 0.38 0.67
t=3 -4.99 -0.25 1.33 0.83
t=4 -4.03 -0.36 1.25 0.61
Lokaler Likeliho o dansatz
0
t
T reat
t
B ase
t
Ag e
t
t=1 -4.03 (0.55) -0.32 (0.06) 1.18 (0.04) 0.64 (0.14)
t=2 -0.79 (0.94) -0.25 (0.09) 0.91 (0.06) -0.05 (0.28) 0.61 (0.04)
t=3 -2.94 (0.55) -0.33 (0.05) 1.11 (0.06) 0.38 (0.13)
t=4 -2.04 (0.69) -0.43 (0.10) 1.03 (0.06) 0.17 (0.19)
79
5.2 Verweildauermo delle
5.2.1 Dauer der Arb eitslosigkeit
Datensatzb eschreibung
Bei diesem Anwendungsb eispiel werden Daten des SOEP (Sozio-
okonomisches
Panel (Hanefeld, 1987)) analysiert. Betrachtet werden insgesamt 1.188 Per-
sonen mit deutscher Staatsangeh
origkeit aus den alten Bundesl
andern, die
in den Jahren 1983-1992 arb eitslos gemeldet waren. Als Einugr
oen wird
das Alter der Person zu Beginn der Arb eitslosigkeit sowie das Geschlecht
aufgenommen.
Mo dellansatz
Die Dauer der Arb eitslosigkeit (in Monaten) wird durch ein Verweildauer-
mo dell mit der Hazardrate
(
t
j
x
it
)=
P
(
T
=
t
j
T
t
;
x
it
)=
exp(
it
)
1+exp(
it
)
und dem linearen Pr
adiktor
it
=
x
0
it
(
t
)+
b
i
mo delliert. Der eindimensionale zuf
allige Eekt
b
i
mo delliert den p erso-
nensp ezischen Einu. Der Regressorenvektor
x
einh
alt neb en dem Inter-
zept als weitere Einugr
oen das Alter und das Geschlecht, d. h. es ist
x
0
it
= (1
;x
SEX
i
;x
AGE
i
). Das Geschlecht ist bin
ar kodiert mit
x
SEX
= 1 f
ur
Frauen und
x
SEX
=0 f
ur M
anner. Die metrische Gr
oe des Alter zu Beginn
der Arb eitslosigkeit geht mit
x
AGE
=(
AGE
35) in das Mo dell ein.
Wir verwenden den Mo dellansatz f
ur Regressionssplines mit kubisch-linearem
Verlauf. Durch den kubischen Verlauf kann die Variabilit
at am Anfang der
Arb eitslosigkeit sehr gut mo delliert werden. Da insb esondere f
ur
t>
20 nur
80
no chsehr wenige Beobachtungen vorliegen, k
onnen durch die Verwendung
des linearen Funktionsansatzes Sch
atzprobleme vermieden werden.
Ergebnisse
In der folgenden Abbildung 5.5 werden die Ergebnisse des kubisch-linearen
Regressionsansatzes mit der Knotenwahl
=18 dem Lokalen Likeliho o dan-
satzes gegen
ub ergestellt.
In b eiden Ans
atzen werden ungef
ahr die gleichen zeitabh
angigen Tendenzen
aufgezeigt, solange die Arb eitslosigkeit nicht mehr als 18 Monate andauert.
Da nach diesem Zeitpunkt die Datenbasis sehr d
unn ist, sollen die Parame-
tersch
atzung f
ur
t>
20 nicht interpretiert werden.
Es wird jedo chauch deutlich, welchen Vorteil der exible Ansatz der Lokalen
Likeliho o d im Gegensatz zur kubisch-linearen Regressions bietet. Zum Bei-
spiel kann der zeitabh
angige Verlauf des Ko eÆzienten
SEX
nicht ad
aquat
durch eine kubisch-lineare Funktion ausgedr
uckt werden.
Das Minimum des kubischen Verlaufes zum Zeitpunkt
t
= 6 mit dem Pa-
rametersch
atzwert
2, entspricht nicht dem zeitabh
angigen Verhalten des
Ko eÆzienten
SEX
. Denn b etrachtet man den Verlauf dieses Ko eÆzien-
ten b eim Lokalen Likeliho o d Ansatz, so werden hier zwei lokale Minima b ei
t
=6 und
t
=14 und vor allem ein globales Minimum f
ur
t
=24 aufgezeigt.
Diese zeitabh
angige Struktur kann nicht durch einen kubisch-linearen Ver-
lauf dargestellt werden, so da das Minimum des Ko eÆzienten
SEX
beim
kubisch-linearen Regressionsansatz nicht durch die Datenstruktur, sondern
durch durch den parametrischen Mo dellansatz b edingt wird.
81
Months
5 101520
-5 -4 -3 -2
INTERCEPT
Months
5 101520
-2 -1 0 1
SEX
Months
5 101520
-0.15 -0.10 -0.05 0.0
AGE
Months
5 101520
-5 -4 -3 -2
INTERCEPT
Months
5 101520
-2 -1 0 1
SEX
Months
5 101520
-0.15 -0.10 -0.05 0.0
AGE
Abbildung 5.5: Parametersch
atzungen mit kubisch-linearem Regressionsan-
satz (links) und Lokalem Likeliho o dansatz (rechts)
5.2.2 Magenkrebsstudie
Datensatzb eschreibung
Der Datensatz eine Magenkrebsstudie (Fahrmeir, 1994) enth
alt die
Ub er-
leb enszeiten (in Monaten) von 90 Patienten einer Magenkrebsstudie. Die
Patienten wurden mit einer Chemotherapie bzw. mit einer kombinierten
Therapieform b ehandelt. Durch die Studie soll der Einu der Therapiefor-
men im Zeitverlauf untersucht werden.
82
Mo dellansatz
Wie b ereits im vorherigen Beispiel b etrachten wir ein Verweildauermo dell
mit der Hazardrate
(
t
j
x
it
)=
P
(
T
=
t
j
T
t
;
x
it
)=
exp(
it
)
1 + exp(
it
)
und dem linearen Pr
adiktor
it
=
0
(
t
)+
x
0
it
TREATMENT
(
t
)+
b
i
:
Der p ersonensp ezischen Einu wird durchden eindimensionale zuf
allige
Eekt
b
i
ber
ucksichtigt. Neb en der Therapieform werden keine weiteren
Einugr
oen b etrachtet. Der Einu der Therapie ist mit
x
= 0 f
ur die
Chemotherapie und mit
x
=1 f
ur die kombinierte Therapieform co diert.
Da in diesem Beispiel die
Ub erleb ensraten der Patienten untersucht werden,
liegen gegen Ende des Beobachtungszeitraumes nur no ch wenige Datens
atze
vor. Wir b etrachten im folgenden Mo delle, die diese Struktur b er
ucksichtigen.
Beim Lokalen Likeliho o d Ansatz verwenden wir den Gaukern und w
ahlen
den Gl
attungsparameter in Abh
angigkeit der Anzahl, der unter Risikoste-
henden Personen, mit
t
=
1
p
n
t
:
Dadurch wird gegen Ende des Beobachtungszeitraumes, in dem nur no ch
wenige Beobachtungen vorhanden sind, eine gr
oere Fensterbreite gew
ahlt
als zu Beginn.
Wir vergleichen die Ergebnisse mit dem Mo dellansatz durch kubisch-lineare
Regression. Dab ei b etrachten wir sowohl eine kubische Funktion, die linear
fortgesetzt wird, als auchden Fall wenn zwei o der drei st
uckweise kubische
Funktionen linear fortgesetzt werden.
83
In den folgenden Abbildungen sind die Hazardraten mit den B
andern der
punktweisen Standardabweichungen f
ur die b eiden Behandlungsarten einge-
zeichnet. Die Hazardraten der Patienten, die mit der kombinierten Therapie-
form b ehandelt wurden (durchgezogene Linie), sind zun
achst h
oher, jedo ch
nach 2 Jahren ist ein gegenteiliger Eekt zu beobachten.
Beim Ansatz der Regressionssplines wurden ein, zwei und drei st
uckweise
kubische Funktionen mit einer linearen Funktion fortgesetzt. Die Stellen der
Knoten sind durch senkrechte Linien gekennzeichnet.
Time (months)
Hazard function
0 1020304050
0.0 0.05 0.10 0.15
gamma: 20
no randomeffect
Time (months)
Hazard function
0 1020304050
0.0 0.05 0.10 0.15
gamma: 30
no randomeffect
Time (months)
Hazard function
0 1020304050
0.0 0.05 0.10 0.15
gamma: 20
with randomeffect
Time (months)
Hazard function
0 1020304050
0.0 0.05 0.10 0.15
gamma: 30
with randomeffect
Abbildung 5.6: Hazardraten b eim Lokalen Likeliho o d Ansatz
84
piecewise-cubic linear regresion splines
no randomeffects
Time (months)
Hazaed function
0 1020304050
0.0 0.05 0.10 0.15
one node
Time (months)
Hazard function
0 1020304050
0.0 0.05 0.10 0.15
two nodes
Time (months)
Hazard function
0 1020304050
0.0 0.05 0.10 0.15
three nodes
Abbildung 5.7: Hazardraten b eim Regressionsspline ohne Ber
ucksichtigung
von zuf
alligen Eekten
85
piecewise-cubic linear regresion splines
with randomeffects
Time (months)
Hazaed function
0 1020304050
0.0 0.05 0.10 0.15
one node
Time (months)
Hazaed function
0 1020304050
0.0 0.05 0.10 0.15
two nodes
Time (months)
Hazaed function
0 1020304050
0.0 0.05 0.10 0.15
three nodes
Abbildung 5.8: Hazardraten b eim Regressionsspline mit Ber
ucksichtigung
von zuf
alligen Eekten
86
Kapitel 6
Anhang
6.1 EM-Algorithmus mit direkter
Gau-Hermite Integrationstechnik
In diesem Abschnitt geb en wir einige technische Details f
ur das verwendete
Verfahren des EM-Algorithmus mit direkter Gau-Hermite Integrationstech-
nik an.
6.1.1 Numerische Integration mit Gau-Hermite
Ist ein Integral der Form
Z
h
(
x
)
dx;
mit
h
(
x
)=
f
(
x
)
g
(
x
)
nicht analytischl
osbar, so kann es durchnumerische Integration (Crourch&
Spiegelman, 1990) approximiert werden. Wir b eschr
anken uns im folgenden
zun
achst auf die Darstellung eines eindimensionalen Integrales, um die No-
tation
ub ersichtlicher zu gestalten.
87
Bei der numerischen Integration wird das Integral durch eine gewichtete Sum-
me approximiert:
Z
R
h
(
x
)
dx
=
Z
R
f
(
x
)
g
(
x
)
dx
M
X
j
=1
f
(
x
j
)
g
j
:
Dab ei bezeichnen
g
j
die Gewichte,
x
j
die St
utzstellen sowie
M
die Anzahl
der St
utzstellen f
ur die numerische Integration.
In unserem Fall entspricht
g
(
x
) der Standardnormalverteilung, d. h. die
betrachtete Funktion hat die Form
h
(
x
)=
f
(
x
)
1
p
2
exp
x
2
2
:
F
ur die Anwendung der Gau-Hermite Integration (Abramowitz & Stegun,
1972), die Funktionen der Form
h
gh
(
x
)=
f
(
x
)
exp
f
x
2
g
b ehandelt, m
ussen die Gewichte und St
utzstellen transformiert werden. Des-
halb substituieren wir
z
=
x
p
2
bzw.
dx
=
p
2
dz
und erhalten
Z
R
h
(
x
)
dx
=
Z
R
f
(
x
)
1
p
2
exp
x
2
2
dx
=
Z
R
1
p
2
p
2
exp
(
z
2
)
f
(
p
2
z
)
dz
M
X
j
=1
1
p
g
j
f
(
p
2
x
j
)
:
Die Gewichte
g
j
und die St
utzstellen
x
j
der Gau-Hermite Integration, welche
tab ellarisiert vorliegen, werden mit
v
j
=
1
p
g
j
und
d
j
=
p
2
x
j
transformiert,
88
so da das Integral mit
Z
R
h
(
x
)
dx
M
X
j
=1
v
j
f
(
d
j
)
approximiert wird.
Die
Ub ertragung auf den mehrdimensionalen Fall erfolgt analog und wir er-
halten
Z
R
q
h
(
x
)
dx
M
X
j
=1
v
j
f
(
d
j
) mit
x
=(
x
1
;:::;x
q
)
0
Dab ei ist
j
2 f
j
1
;:::;j
q
g
ein Multiindex. Es wird insgesamt
ub er
M
=
M
1
:::
M
q
St
utzstellen approximiert. Die Gewichte
v
j
= (
v
j
1
;:::;v
j
q
)
und die zugeh
origen St
utzstellen
d
j
= (
d
j
1
;:::;d
j
q
)sind Elemente aus den
Kartesischen Pro dukten
v
j
2
V
1
V
q
bzw.
d
j
2
D
1
D
q
:
Dab ei ist
V
l
=
f
v
1
;:::;v
M
l
g
die Menge der transformierten Gewichte und
D
l
=
f
d
1
;:::;d
M
l
g
die Menge der transformierten St
utzstellen bei der Wahl
von
M
l
St
utzstellen f
ur
l
=1
;:::;q
.
Das b edeutet, da die Anzahl der St
utzstellen bei mehrdimensionalen Inte-
gralen sehr gro wird. Soll z. B. ein dreidimensionales Integral approximiert
werden, wob ei in jeder Dimension 10 St
utzstellen verwendet werden, so wer-
den insgesamt 10
10
10 = 1000 St
utzstellen verwendet.
6.1.2 Sch
atzung der festen Parameter
Die Sch
atzung der festen Parameter erfolgt durch den EM-Algorithmus unter
Verwendung der ob en b eschrieb enen Gau-Hermite Integrationstechnik.
89
Dichte
Die Beobachtungseinheiten
y
i
=(
y
i
1
;:::;y
in
i
)sindf
ur
i
=1
;:::;N
,gegeben
,
Q
und
b
1
;:::;b
N
, b edingt unabh
angig. Damit gilt f
ur die Dichte von
Y
=(
y
1
;:::;y
N
)
0
f
(
Y
j
b
1
;:::;b
N
;
)=
N
Y
i
=1
f
(
y
i
j
b
i
;
)mit
f
(
y
i
j
b
i
;
)=
n
i
Y
t
=1
f
(
y
it
j
b
i
;
)
:
Die Sch
atzung der Parameter
und
Q
kann durch Maximierung der Margi-
nalen Log-Likeliho o d
l
(
; Q
)=
N
X
i
=1
ln
L
i
(
; Q
)mit
L
i
(
; Q
)=
Z
f
(
y
i
j
b
i
;
)
p
(
b
i
;Q
)
db
i
durchgef
uhrt werden. Dab ei b ezeichnet
p
(
b
i
;Q
) die Dichte von
b
i
,die als
normalverteilt vorausgesetzt wurde.
Das Integral
ub er diese Mischungsdichte kann nur f
ur einige Sp ezialf
alle ana-
lytisch gel
ost werden, so da f
ur die L
osung des Integrales numerische Ver-
fahren verwendet werden.
Bevor diese Verfahren angewendet werden k
onnen, wird der Mo dellansatz
f
ur den linearen Pr
adiktor
it
=
Z
it
+
W
it
b
i
umparametrisiert.
Umparametrisierung des Linearen Pr
adiktors
F
ur die zuf
alligen Eekte
b
i
hatten wir vorausgesetzt, da
b
i
iid
N
(0
;Q
) gilt.
Da
Q
eine Kovarianzmatrix ist, also p ositiv denit ist, l
at sich
Q
durch
die Cholesky-Zerlegung in
Q
=
Q
1
=
2
Q
T=
2
zerlegen, wob ei
Q
1
=
2
eine untere
90
Dreiecksmatrix ist und
Q
T=
2
die dazu transp onierte ob ere Dreiecksmatrix.
Die ob ere Dreiecksmatrix
Q
1
=
2
l
at sich vektorisieren mit:
=(
q
11
;:::;q
1
r
;q
22
;:::;q
2
r
;:::;q
rr
)
0
mit
0
B
B
@
q
11
q
12
::: q
1
r
0
q
22
::: q
2
r
.
.
.
:::
.
.
.
q
rr
1
C
C
A
=
Q
1
=
2
:
Der Vektor
=
vec
(
Q
1
=
2
)enth
alt nur die Ko eÆzienten aus
Q
1
=
2
, welche
nicht nach Konstruktion von
Q
1
=
2
gleich 0 sind. Ist
bekannt, so ist auch
die Kovarianzmatrix
Q
eindeutig b estimmt.
Mit der Umparametrisierung
b
i
=
Q
1
=
2
a
i
und
a
i
iid
N
(0
;
I)
;
erhalten wir
W
it
b
i
=
W
it
Q
1
=
2
a
i
=[
a
0
i
W
it
]
unter Verwendung des Kroneckerpro duktes
.
Damit folgt f
ur den linearen Pr
adiktor
it
= [
Z
it
; a
0
i
W
it
]
mit
=
!
und
a
i
iid
N
(0
;
I).
Diesen Mo dellansatz f
ur den linearen Pr
adiktor werden wir im folgenden
verwenden, denn er hat zwei entscheidende Vorteile. Zum einen wurde der
lineare Pr
adiktor in die
ubliche Form eines linearen Mo dells
ub erf
uhrt. Durch
eine Sch
atzung von
erh
alt man Sch
atzer f
ur die b eiden festen Parameter
und
Q
. Der zweite Vorteil ist, da in der Mischungsdichte von
Y
die Dichte
p
(
b
i
;Q
), welche eine normalverteilte Dichte mit unbekannter Kovarianz ist,
durch eine standardnormalverteilte Dichte, d. h. durch
p
(
a
i
;
I) ersetzt wer-
den kann.
91
Sch
atzung von
und
Q
Das Mo dell
it
=
Z
it
+
W
it
b
i
mit
b
i
N
(0
;Q
) liegt nach der Parametrisie-
rung in der Form:
it
=[
Z
it
; a
0
i
W
it
]
mit
a
i
iid
N
(0
;
I)
mit
=
!
und
=
vec
(
Q
1
=
2
)mit
Q
=
Q
1
=
2
Q
T=
2
vor.
Falls die
a
i
bekannt sind, ist das Problem ein gew
ohnliches Generalisier-
tes Lineares Mo dell. Da dies ab er nicht der Fall ist, wenden wir den EM-
Algorithmus an.
Der EM-Algorithmus, der auch Grundlage des EM-Typ e Algorithmus ist,
derinAbschnitt 2.1.2 auf Seite 25 vorstellt wurde ist ein iteratives Verfahren
(Dempster, Laird & Rubin, 1977; Little & Rubin, 1987; McLachlan & Kris-
hnan, 1997; Sundb erg, 1974) f
ur die Bestimmung eines Maximum Likeliho o d
Sch
atzers in unvollst
andigen Datensituationen. Der Algorithmus b esteht in
jedem Iterationsschritt aus einem
Expectation-Schritt
(E-Schritt) und einem
Maximizing-Schritt
(M-Schritt).
E-Schritt:
(Berechnung des b edingten Erwartungswertes)
M
(
j
h
k
i
)=
E
ln
f
(
Y; A
j
)
j
Y;
h
k
i
mit
Y
=
0
B
B
@
y
1
.
.
.
y
N
1
C
C
A
und
A
=
0
B
B
@
a
1
.
.
.
a
N
1
C
C
A
.
92
Anschlieend wird diese Log-Likeliho o d im M-Schritt maximiert, um einen
neuen Sch
atzer
h
p
+1
i
f
ur
zu erhalten.
M-Schritt:
(Maximierung von
M
(
;
h
k
i
))
@M
(
;
h
k
i
)
@
=0
Im EM-Algorithmus wird im E-Schritt der b edingte Erwartungswert
M
(
j
h
k
i
) =
E
ln
f
(
Y; A
j
)
j
Y;
h
k
i
mit
Y
=
0
B
B
@
y
1
.
.
.
y
N
1
C
C
A
und
A
=
0
B
B
@
a
1
.
.
.
a
N
1
C
C
A
b erechnet.
Da wir im M-Schritt diesen Erwartungswert bez
uglich
maximieren wollen,
sind f
ur uns konstante Terme o der Terme, die von
unabh
angig sind, ohne
Bedeutung.
Zun
achst ist
M
(
j
h
k
i
) =
E
ln
f
(
Y; A
j
)
j
Y;
h
k
i
=
Z
ln
f
(
Y; A
j
)
f
(
A
j
Y;
h
k
i
)
dA :
Der erste Term des obigen Integrals lautet
93
ln
f
(
Y; A
j
) = ln (
f
(
Y
j
A;
)
g
(
A
))
=
N
X
i
=1
ln
f
(
y
i
j
a
i
;
)+
N
X
i
=1
ln
g
(
a
i
)
;
wob ei
N
P
i
=1
ln
g
(
a
i
) unabh
angig von
ist.
F
ur den zweiten Term ist
f
(
A
j
Y;
h
k
i
) =
f
(
Y; A
j
h
k
i
)
f
(
Y
j
h
k
i
)
=
f
(
Y
j
A;
h
k
i
)
g
(
A
)
R
f
(
Y
j
A;
h
k
i
)
g
(
A
)
dA
=
N
Q
i
=1
f
(
y
i
j
a
i
;
h
k
i
)
N
Q
i
=1
g
(
a
i
)
N
Q
i
=1
R
f
(
y
i
j
a
i
;
h
k
i
)
g
(
a
i
)
da
i
:
Der Nenner, der von
abh
angt, ist eine Konstante, die f
ur die Maximierung
nicht ber
ucksichtigt werden mu.
Die zu maximierende Funktion l
at sich nun vereinfachen zu
~
M
(
j
h
k
i
) =
N
X
i
=1
Z
ln
f
(
y
i
j
a
i
;
)
f
(
y
i
j
a
i
;
h
k
i
)
g
(
a
i
)
da
i
:
Durch Gau-Hermite Integration wird das Integral durch folgende Summe
approximiert
~
M
(
;
h
k
i
)
M
GH
(
;
h
k
i
)=
N
X
i
=1
M
X
j
=1
ln
f
(
y
i
j
d
j
;
)
v
j
f
(
y
i
j
d
j
;
h
k
i
)
M
P
s
=1
v
j
f
(
y
i
j
d
s
;
h
k
i
)
94
Die Gewichte
v
j
und die St
utzstellen
d
j
sind die transformierten Gewichte
bzw. St
utzstellen der Gau-Hermite Integration, wie sie im vorherigen Ab-
schnitt deniert wurden.
Wir denieren
c
ij
=
v
j
f
(
y
i
j
d
j
;
h
k
i
)
M
P
s
=1
v
s
f
(
y
i
j
d
s
;
h
k
i
)
und maximieren im M-Schritt
M
GH
(
j
h
k
i
), durch Nullsetzen der Ableitung
M
GH
(
j
h
k
i
), d. h.
@M
GH
(
j
h
k
i
)
@
=
N
X
i
=1
M
X
j
=1
c
ij
@
ln
f
(
y
i
j
d
j
;
)
@
!
=0
Da
=(
;
)
0
ist, folgt
@
ln
f
(
y
i
j
d
j
;
)
@
=
0
B
B
@
@
ln
f
(
y
i
j
d
j
;
)
@
@
ln
f
(
y
i
j
d
j
;
)
@
1
C
C
A
und damit entspricht
@M
GH
(
j
h
k
i
)
@
der Scorefunktion dem mit
c
ij
gewichtete GLM:
itj
=
Z
it
; d
0
j
W
it
!
mit
h
(
itj
)=
E
(
y
it
) f
ur
j
=1
;:::;M
.
Die Parameter
und
werden mit einem gewichteten GLM gesch
atzt, wo-
bei die Designmatrix
Z
it
; d
0
j
W
it
neb en den Regressoren auch no chdie
St
utzstellen f
ur die Gau-Hermite Integration enth
alt.
Die Resp onsematrix wird mit
y
itj
=
y
it
f
ur
j
= 1
;:::;M
entsprechend an-
gepat. Jede Beobachtung
y
itj
wird f
ur
t
= 1
;:::;n
i
bei der Sch
atzung des
95
GLM's mit
c
ij
=
c
ij
(
h
k
i
)=
v
j
f
(
y
i
j
d
j
;
h
k
i
)
M
P
s
=1
v
s
f
(
y
i
j
d
s
;
h
k
i
)
gewichtet.
96
6.2 EM-Typ e Algorithmus
6.2.1 Blo ckdiagonalgestalt der Fishermatrix
Da b eim EM-Typ e Algorithmus die Fishermatrix in Blo ckdiagonalgestalt
vorliegt, kann diese Struktur f
ur die Parametersch
atzung mit dem Fisher-
Scoring Algorithmus ausgenutzt werden.
Mit
F
=
2
6
6
4
F
1
1
O
.
.
.
O
F
T
T
3
7
7
5
und
F
b
i
=
2
6
6
4
F
1
b
i
.
.
.
F
T
b
i
3
7
7
5
=
F
0
b
i
folgt:
F
(
Æ
)
h
p
i
=
2
6
6
6
6
6
4
F
F
b
1
F
b
N
F
b
1
F
b
1
b
1
O
.
.
.
.
.
.
F
b
N
O
F
b
n
b
n
3
7
7
7
7
7
5
Ist
Æ
h
p
+1
i
=
Æ
h
p
i
+
4
Æ
h
p
i
;
so folgt b eim Fisher-Scoring Algorithmus:
Æ
h
p
+1
i
=
Æ
h
p
i
+
F
1
Æ
h
p
i
s
Æ
h
p
i
;
da
s
Æ
h
p
i
=
F
Æ
h
p
i
4
Æ
h
p
i
gilt.
Da die Fishermatrix Blo ckdiagonalgestalt hat, gilt
s
h
p
i
=
F
4
h
p
i
+
N
X
i
=1
F
b
i
4
b
h
p
i
i
s
b
h
p
i
i
=
F
b
i
4
h
p
i
+
F
b
i
b
i
4
b
h
p
i
i
97
Au
osen der 2. Gleichung nach
4
b
h
p
i
i
f
uhrt zu
4
b
h
p
i
i
=
F
1
b
i
b
i
n
s
b
h
p
i
i
F
b
i
4
h
p
i
o
und mit Einsetzen von
4
b
h
p
i
i
in die 1. Gleichung folgt
s
h
p
i
=
F
4
h
p
i
+
N
X
i
=1
F
b
i
F
1
b
i
b
i
n
s
b
h
p
i
i
F
b
i
4
h
p
i
o
=
F
4
h
p
i
+
N
X
i
=1
F
b
i
F
1
b
i
b
i
s
b
h
p
i
i
4
h
p
i
N
X
i
=1
F
b
i
F
1
b
i
b
i
F
b
i
:
So da f
ur
4
h
p
i
gilt:
4
h
p
i
=
(
F
N
X
i
=1
F
b
i
F
1
b
i
b
i
F
b
i
)
1
(
s
h
p
i
N
X
i
=1
F
b
i
F
1
b
i
b
i
s
b
h
p
i
i
)
:
98
Zusammenfassung
In vielen Studien werden die Beobachtungseinheiten
ub er einen l
angeren Zeit-
raum b etrachtet. Es liegt somit einerseits eine grupp en- und/o der individu-
ensp ezische Datenstruktur vor und andererseits mu auch ein m
oglicher
zeitlicher Einu b er
ucksichtigt werden. In dieser Arb eit wurden verschiede
Mo dellans
atze vorgestellt, die den grupp en- bzw. individuensp ezischen Ein-
u durch die Aufnahme von zuf
alligen Eekten b er
ucksichtigen und gleich-
zeitig einen zeitlichen Einu durch die Betrachtung von zeitvariierenden
Ko eÆzienten erm
oglichen.
In den Simulationsstudien und Anwendungsb eispielen hat sich gezeigt, da
der Lokale Likeliho o dansatz, der durch die Aufnahme von zuf
alligen Eekten
erweitert wurde, von allen vorgestellten Mo dellans
atzen am besten geeignet
ist.
Bei diesem Ansatz werden durch eine Kernfunktion die einzelnen Beobach-
tungen in Abh
angigkeit des b etrachteten Beobachtungszeitpunktes unter-
schiedlich gewichtet, so da durch die Wahl des Gl
attungsparameters der
Kernfunktion ein sehr exibler Gl
attungsansatz m
oglich ist. Der optima-
le Gl
attungsparameter kann durch Kreuzvalidierung, die jedo ch mit groem
Rechenaufwand verbunden ist, b estimmtwerden.
Selbstverst
andlichm
ussen die Varianzen der Parametersch
atzungen f
ur deren
Interpretation immer ber
ucksichtigt werden. Jedo ch k
onnen durch die Be-
99
trachtung der punktweisen Kondenzintervalle alle Parametersch
atzungen,
die in Abh
angigkeit des Gl
attungsparameters erzielt wurden, inhaltlich in-
terpretiert werden, so da auch Parametersch
atzungen interpretiert werden
k
onnen, die nicht mit dem optimalen Gl
attungsparameter b estimmt wur-
den. Damit ist es durch den Lokalen Likeliho o dansatz m
oglich, einerseits
durch eine geringe Gl
attung die einzelnen betrachteten Zeitpunkte separat
zu analysieren und andererseits durcheine starke Gl
attung einen globalen
Trend
ub erblick
ub er den gesamten Zeitraum zu erhalten.
Beispiele f
ur die Auswirkungen unterschiedlicher Gl
attungsparameter wur-
den im Emesis Datensatz in Abschnitt 5.1.1 gezeigt.
Beim Lokalen Likeliho o dansatz mit zuf
alligen Eekten wurde im Abschnitt
3.1.2 zwei Verfahren, der EM-Algorithmus mit Gau-Hermite Integrations-
technik und der EM-Typ e Algorithmus, f
ur die Bestimmung der Parame-
tersch
atzwerte angegeb en. Der EM-Algorithmus mit Gau-Hermite Inte-
grationstechnik sollte trotz seines groen Rechenaufwands aufgrund der be-
trachteten Integrationsst
utzstellen, bevorzugt angewendet werden, da bei
diesem Verfahren auchbei einer deutlichen Auspr
agung der Varianz der
grupp en- bzw. individuensp ezischen Ein
usse die Parametersch
atzungen
durchgef
uhrt werden k
onnen. Im Simulationsb eispiel 4.1.3 wurde gezeigt,
da der EM-Typ e Algorithmus in diesem Fall nicht geeignet ist, da die
Sch
atzungen der zuf
alligen Eekte unzureichend sind.
F
ur den EM-Algorithmus mit Gau-Hermite Integrationstechnik, konnte im
Simulationsb eispiel 4.1.2 als auch in den Anwendungsb eispielen in Kapi-
tel 5 gezeigt werden, da dieses Verfahren gut geeignet ist, die Parame-
tersch
atzungen selbst und die Bestimmung ihrer Varianzen durchzuf
uhren.
F
ur den in Abschnitt 2.1.3 vorgestellten Mo dellansatz der Nicht Parame-
trischen Maximum Likeliho o d konnte keine zufriedenstellende Erweiterung
100
durch die Ber
ucksichtigung von zeitvariierenden Ko eÆzienten entwickelt wer-
den.
Ein weiterer sehr einfacher Ansatz f
ur die Ber
ucksichtigung von zeitvariieren-
den Ko eÆzienten und zuf
alligen Eekten ist der Marginale Likeliho o dansatz.
Bei diesem Ansatz kann durchden mo dellb edingten Aufbau der Designma-
trix nur einen geringe Anzahl von b etrachteten Zeitpunkten ber
ucksichtigt
werden, da sonst in der Regel numerische Probleme b ei der Bestimmung der
Parametersch
atzungen entstehen. Deshalb kann dieser Mo dellansatz in der
Praxis nur sehr b eschr
ankt eingesetzt werden.
Der Vorteil des exiblen Lokalen Likeliho o dansatzes gegen
ub er der parame-
trischen und semi-parametrischen Ans
atze wurde in den Anwendungsb ei-
spielen f
ur die Dauer der Arb eitslosigkeit in Abschnitt 5.2.1 und der Ma-
genkrebsstudie in Abschnitt 5.2.2 dargestellt werden. Wird keine dem Mo-
dell entsprechende parametriche Form gew
ahlt,sok
onnen sehr leicht falsche
R
uckschl
usse getroen werden.
Zus
atzlich wurde in Abschnitt 5.1.2 am Beispiel der epileptischen Anf
alle
deutlich, da die Bestimmung der Varianzen der Parametersch
atzungen mit
der Delta-Metho de b ei der Polynomialen Regression aus numerischen Gr
unden
nicht immer m
oglich ist.
101
Literatur
Abramowitz, M. und Stegun, I. (1972).
Handbook of Mathematical Func-
tions
.New York: Dover.
Agresti, A. (1990).
Categorical Data Analysis
.NewYork: Wiley.
Aitkin, M. (1995). NPML estimation of the mixing distribution in general
statistical mo dels with unobserved random eects. (
working-paper
).
Aitkin, M. (1996). A general maximum likeliho o d analysis of overdisp ersi-
on in generalized linear mo dels.
Statistics and Computing 6
, 251{262.
Aitkin, M. und Aitkin, I. (1996). A hybrid EM/Gauss-Newton algorithm
for maximum likeliho o d in mixture distributions. (
working-paper
).
Aitkin, M. und Francis, B. (1995). Fitting overdisp ersed generalized linear
mo dels by nonparametric maximum likeliho o d.
GLIM Newsletters 25
,
37{45.
Aitkin, M. und Wilson, G. (1980). Mixture mo dels, outliers, and the EM
algorithm.
Technometrics 22
(3), 325{331.
Bo oth, J. und Hob ert, J. (1999). Maximizing generalized linear mixed
mo dels with an automated monte carlo em alorithm.
Journal of the
Royal Statistical Society Ser. B 61
, 265{285.
Breslow, N. und Clayton, D. (1993). Approximate inference in generalized
linear mixed mo dels.
Journal of the American Statistical Associati-
102
on 88
(421), 9{25.
Crourch, E. und Spiegelman, D. (1990). The evaluation of integrals of
the form
R
+
1
1
f
(
t
)
exp
(
t
2
)
dt
: Application to logistic-normal mo dels.
Journal of the American Statistical Association 85
(410), 464{469.
Davis, C. (1991). Semi-parametric and non-parametric metho ds for the
analysis of rep eated measurements with applications to clinical trials.
Statistics in Medicine 10
, 1959{1980.
Dempster, A., Laird, N., und Rubin, D. (1977). Maximum likeliho o d from
incomplete data via the EM-algorithm.
Journal of the Royal Statistical
Society B 39
,1{38.
Dietz, E. und B
ohning, D. (1994). Analysis of longitudinal data using nite
mixture mo dels.
Statistische Hefte 35
, 203{210.
Diggle, P., Liang, K., und Zeger, S. (1994).
Analysis of Longitudinal Data
.
Oxfort University Press.
Efron, B. (1988). Logistic regression, survival analysis, and the Kaplan-
Meier-curve.
Journal of the American Statistical Association 83
,414{
425.
Eilers, P. und Marx, B. (1996). Flexible smo othing with B-splines an p en-
alties.
Statistical Science 11
(2), 89{121.
Eubank, R. (1988).
Spline Smoothing and Nonparametric Regression
.New
York: Marcel Dekker.
Everitt, B. und Hand, D. (1981).
Finite Mixture Distributions
. London:
Chapman and Hall.
Fahrmeir, L. (1994). Dynamic mo delling and p enalized likeliho o d estima-
tion for discrete time survival data.
Biometrika 81
(2).
Fahrmeir, L., Hamerle, A., und Tutz, G. (1994).
Multivariate statistische
Verfahren
. Berlin: de Gruyer.
103
Fahrmeir, L. und Tutz, G. (1994).
Multivariate Statistical Model ling Based
on Generalized Linear Models
. New York: Springer Verlag.
Fan, J., Hall, P., Martin, M., und Patil, P. (1996). On lo cal smo othing
of nonparametric curve estimators.
Journal of the American Statistical
Association 91
(433).
Friedman, J. (1991). Multivariate adaptive regression splines.
The Annals
of Statistics 19
, 1{141.
Geyer, C. und Thompson, E. (1992). Constrained monte carlo maximum
likeliho o d for dep end data.
Journal of the Royal Statistical Society Ser.
B 54
, 657{683.
Gourieroux, C. und Montfort, A. (1989). Simulation based infereerce in
mo dels with heterogeneity.
Document de Trvail INSEE/ENSAE No.
8902
.
Hanefeld, U. (1987).
Das sozio-
okonomische Panel
.Frankfurt: Campus.
Hastie, T. und Loader, C. (1993). Lo cal regression: Automatic kernel car-
pentry.
Statistical Science 8
(2), 120{143.
Hastie, T. und Tibshirani, R. (1990).
GeneralizedAdditive Models
. London:
Chapman and Hall.
Hastie, T. und Tibshirani, R. (1993). Varying-co eÆcient mo dels.
Journal
of the Royal Statistical Society B 55
(4), 757{796.
Hennevogl, W. (1991).
Sch
atzung generalisierter Regressions- und Zeitrei-
henmodel le mit varriierenden Parametern
. Dissertation, Universit
at
Regensburg.
Hinde, J. (1982). Comp ound p oisson regression mo dels. In R. Gilchrist
(Hrsg.),
GLIM'82 Internat. Conf. Generalized Linear Models
,New
York, S. 109{121. Springer.
104
Jansen, R. (1993). Maximum likeliho o d in a generalized linear nite mix-
ture mo del by using th EM algorithm.
Biometrics 49
, 227{231.
Kauermann, G. und Tutz, G. (1995). Lo cal likeliho o d estimation and bias
reduction in varying co eÆcient mo dels. Forschungsb ericht, Technische
Universit
at Berlin. 95-9.
Kauermann, G. und Tutz, G. (2000). Lo cal likeliho o d estimation in
varying-co eÆcient mo dels including additive bias correction.
Journal
of Nonparametric Statistics 12
, 343{371.
Laird, N. M. und Louis, T. A. (1982). Approximate p osterior distributions
for incomplete data problems.
Journal of the Royal Statistical Society B
44
(2), 190{200.
Little, R. und Rubin, D. (1987).
Statistical Analysis with Missing Data
.
New York: John Wily & Sons.
Longford, N. (1993).
Random CoeÆcient Models
. Oxford: Clarendon
Press.
McCullagh, P. (1980). Regression mo del for ordinal data (with discussion).
Journal of the Royal Statistical Society B 42
, 109{127.
McCullagh, P. und Nelder, J. (1989).
Generalized Linear Models
(2 Au.).
London, New York: Chapman and Hall.
McCullo ch, C. (1994). Maximum likeliho o d variance comp onents estima-
tion for binary data.
Journal of the American Statistical Associati-
on 89
(425), 330{335.
McCullo ch, C. (1997). Maximum likeliho o d algorithms for generalized li-
near mixed mo dels.
Journal of the American Statistical Association 92
,
162{170.
McCullo ch, C. (2000). An intro duction to generalized linear mixed mo dels.
(
working-paper
).
105
McLachlan, G. und Basford, K. (1988).
Mixture Models. Inference and
Applications to Clusterung
.New York: Marcel Dekker.
McLachlan, G. und Krishnan, T. (1997).
The EM Algorithm and Extensi-
ons
.New York: Wiley.
Meng, X.-L. (1997). The EM alogithm and medical studies: a historical
link.
Statistical Methods in Medical Research 6
,3{23.
Meng, X.-L. und van Dyk, D. (1997a). The EM algorithm - an old folk-song
sung to a fast new time.
B 59
(3), 511{567.
Meng, X.-L. und van Dyk, D. (1997b). Fast EM-typ e implementations for
mixed-eects mo dels.
Journal of the Royal Statistical Society
. (Submit-
ted).
Nelder, J. A. und Wedderburn, R. W. M. (1972). Generalized linear mo-
dels.
Journal of the Royal Statistical Society A 135
, 370{384.
Quintana, R., Lui, J., und Pino, G. (1999). Monte carlo em with imp ort-
ance reweighting and it application in random eects mo dels.
Compu-
tational Statistics and Data Analysis 29
, 429{444.
Rice, J. und Silverman, B. (1991). Estimating the mean and covariance
structure nonparametrically when the data are curves.
Journal of the
Royal Statistical Society B 53
(1), 233{243.
Rupp ert, D., Reish, R., und Carroll, R. (1984). Optimization using sto-
chastic approximation and monte carlo simulation (with application to
harvesting of atlantic menhaden.
Biometrics 40
, 535{545.
Silverman, B. (1984). Spline smo othing: The equivalent variable kernel
metho d.
The Annals of Statistics 12
(3), 898{916.
Stone, C., Hansen, M., Ko op erb erg, C., und Truong, Y. (1997). Polynomial
splines and their tensor pro ducts in extended linear mo dels.
The Annals
of Statistics 25
, 1371{1470.
106
Sundb erg, R. (1974). Maximum likeliho o d theory for incomplete data from
an exp onential family.
Scand. J. Statist. 1
, 49{58.
Thall, P. und Vail, S. (1990). Some covariance mo dels for longitudinal
count data with overdisp ersion.
Biometrics 46
, 657{671.
Titterington, D., Smith, A., und Makov, U. (1985).
Statistical Analysis of
Finite Mixture Distributions
.New York: Wiley.
Tutz, G. (1989). On cross-validation for discrete kernel estimates in dis-
crimination.
Communications in Statistics, Theory and Methods 11
,
4145{4162.
Tutz, G. (1999). Varying co eÆcients in generalized linear random ef-
fects mo dels: A lo cal likeliho o d approach.
Discussion Paper, Ludwig-
Maximilian-Universit
at M
unchen,Sonderforschungsbereich 386 171
.
Tutz, G. und Hennevogl, W. (1996). Random eects in ordinal regression
mo dels.
Computational Statistics and Data Analysis 22
, 537{557.
Tutz, G. und Kauermann, G. (1995). Varying co eÆcients in multivariate
generalized linear mo dels: a lo cal likeliho o d approach.
Forschungsbe-
reichte des Fachbereiches Informatik, TU-Berlin
.
Tutz, G. und Kauermann, G. (1997). Lo cal estimators in multivariate ge-
neralized linear mo dels with varying co eÆcients.
Computational Stati-
stics 12
, 193{208.
van der Linde, A. (1994). On cross-validation for smo othing splines in the
case of dep endent obsersartions.
36
(1), 67{73.
Wu, C. (1983). On the convergence prop erties of the EM algorithm.
The
Annals of Statistics 11
, 95{103.
107
