Technische Universiteit Delft
Faculteit Elektrotechniek, Wiskunde en Informatica
Delft Institute of Applied Mathematics
Periodieke controle op prostaatkanker: doen of niet doen?
Competing risk analyse bij voorspellingen in de geneeskunde.
(Engelse titel: Periodic screening for prostate cancer:
to do or not to do?)
(Engelse ondertitel: Medical predictions using competing risk analysis.)
Verslag ten behoeve van het
Delft Institute of Applied Mathematics
als onderdeel ter verkrijging
van de graad van
BACHELOR OF SCIENCE
in
TECHNISCHE WISKUNDE
door
D.S. de Groot
Delft, Nederland
Juli 2014
c 2014 door D.S. de Groot. Alle rechten voorbehouden.
Copyright BSc verslag TECHNISCHE WISKUNDE
“Periodieke controle op prostaatkanker: doen of niet doen?”
“Competing risk analyse bij voorspellingen in de geneeskunde.”
(Engelse titel: “Periodic screening for prostate cancer:
to do or not to do?”)
(Engelse ondertitel: “Medical predictions using competing risk analysis.”)
D.S. de Groot
Technische Universiteit Delft
Begeleiders
Dr. A.J. Cabo
Prof.dr.ir. G. Jongbloed
Overige commissieleden
Prof.dr.ir. A.W. Heemink
Dr. B. van den Dries
Juli, 2014
Delft
Samenvatting
Door het ErasmusMC in Rotterdam is een prostaatkanker-dataset beschikbaar gesteld. De
onderzoeksvraag bij deze dataset is ”Is het relevant om te controleren op prostaatkanker bij
mannen tussen de 54 jaar en 75 jaar”. Om deze vraag te beantwoorden is gebruik gemaakt van
verschillende analyses. Bovendien zijn twee datasets gesimuleerd.
Survivalanalyse is een methode om te analyseren hoe lang het duurt tot een gebeurtenis
plaatsvindt. In de geneeskunde is dit een vaak gebruikte analyse. Survivalanalyse wordt vaak
toegepast op data waarvan alleen bekend is dat tot het einde van het onderzoek de gebeurtenis niet plaats gevonden heeft, dit heet rechts-gecensureerde data. Dat de gebeurtenis niet
plaatsvindt tijdens het onderzoek kan meerdere oorzaken hebben. In de geneeskunde komt het
vaak voor dat data rechts-gecensureerd is doordat een persoon overleden is door een andere gebeurtenis, waardoor de onderzochte gebeurtenis niet plaats kan vinden. Dit heet een competing
risk.
Standaard survivalanalyse gaat er vanuit dat het tijdstip van censuur onafhankelijk is van
het tijdstip waarop de gebeurtenis plaatsvindt. Als de gebeurtenis niet plaatsvindt door het
optreden van een competing risk, is deze aanname in het algemeen niet waar. In dit geval moet
gebruik gemaakt worden van competing risk analyse. De twee competing risk modellen die het
meest gebruikt worden zijn het oorzaak-specifieke hazardmodel en het model van Fine en Gray.
Voor de gesimuleerde datasets geldt dat bij dataset 1 de ene groep meer kans heeft op de
beoogde gebeurtenis dan de andere groep. Bij dataset 2 zijn de kansen op gebeurtenissen gelijk.
Na analyse van deze datasets volgt de conclusie dat als er verschil is in het aantal personen
dat overlijdt aan de beoogde gebeurtenis, dit ook naar voren komt uit de resultaten. Voor de
prostaatkanker-data geldt dat de competing risks de doodsoorzaak domineren. Daardoor komt
uit een globale analyse een ander resultaat dan uit een meer nauwkeurige analyse. Uit meer
nauwkeurige analyse volgt dat het relevant is om te controleren op prostaatkanker bij mannen
tussen de 54 jaar en 75 jaar.
5
6
Inhoudsopgave
1 Inleiding
9
2 Survivalanalyse
2.1 Gecensureerde data . . . . . . . . . .
2.2 Functies . . . . . . . . . . . . . . . .
2.3 Methoden . . . . . . . . . . . . . . .
2.3.1 Methoden zonder covariaten .
2.3.2 Methode met covariaten . . .
.
.
.
.
.
11
11
11
13
13
17
3 Competing risk analyse
3.1 Methoden . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.1.1 Methode zonder covariaten . . . . . . . . . . . . . . . . . . . . . . . . . .
3.1.2 Methoden met covariaten . . . . . . . . . . . . . . . . . . . . . . . . . . .
23
24
24
25
4 Analyse van prostaatkanker-data
4.1 Databeschrijving . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.2 Resulaten . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.2.1 Aantal overleden mannen per groep . . . . . . . . . . . . . .
4.2.2 Survivalanalyse . . . . . . . . . . . . . . . . . . . . . . . . . .
4.2.3 Competing risk analyse . . . . . . . . . . . . . . . . . . . . .
4.2.4 Aantal mannen dat overleden is aan prostaatkanker per groep
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
29
29
30
30
31
34
36
5 Simulatie
5.1 Simulatie beschrijving . . . . . .
5.2 Resultaten . . . . . . . . . . . . .
5.2.1 Aantal overleden personen
5.2.2 survivalanalyse . . . . . .
5.2.3 Competing risk analyse .
. . . . . .
. . . . . .
per groep
. . . . . .
. . . . . .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
39
39
39
40
42
48
6 Conclusie
6.1 Conclusies prostaatkanker-data
6.2 Conclusies simulatie . . . . . .
6.2.1 Conclusies dataset 1 . .
6.2.2 Conclusies dataset 2 . .
6.3 Slotconclusie . . . . . . . . . .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
51
51
52
52
53
53
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
7 Discussie
55
A Functies
57
7
8
INHOUDSOPGAVE
B Dataset
63
C Data simuleren
65
D Analyse met verhoudingen
D.1 Prostaatkanker-dataset . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
D.2 Gesimuleerde datasets . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
67
67
68
E Survivalanalyse
E.1 Prostaatkanker-dataset . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
E.2 Gesimuleerde datasets . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
71
71
72
F Competing risk analyse
F.1 Prostaatkanker-dataset . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
F.2 Gesimuleerde datasets . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
75
75
75
Hoofdstuk 1
Inleiding
Periodieke controle op prostaatkanker: doen of niet doen? Om hierop een antwoord te vinden
zullen verschillende analyses aan bod komen.
De survivalanalyse geeft antwoord op de vraag ’hoe lang duurt het voordat. . . ’. Dit is een
veel voorkomende vraag in de geneeskunde. Survivalanalyse wordt vaak toegepast op data waarvan alleen bekend is dat tot het einde van het onderzoek de gebeurtenis niet plaats gevonden
heeft, dit heet rechts-gecensureerde data. Dat de gebeurtenis niet plaatsvindt tijdens het onderzoek kan meerdere oorzaken hebben. In de geneeskunde komt het vaak voor dat de onderzochte
gebeurtenis niet plaatsvindt doordat de persoon overleden is door een andere gebeurtenis. Een
andere gebeurtenis waardoor de onderzochte gebeurtenis niet plaats kan vinden heet een competing risk.
Standaard survivalanalyse gaat er vanuit dat het tijdstip waarop een persoon uit het onderzoek verdwijnt, zonder dat de gebeurtenis plaatsgevonden heeft, onafhankelijk is van het tijdstip
waarop de gebeurtenis plaatsvindt. Als de gebeurtenis niet plaatsvindt door het optreden van
een competing risk is deze aanname in het algemeen niet waar. In dit geval moet gebruik gemaakt worden van competing risk analyse. Als survivalanalyse toch gebruikt wordt, geeft dit
verkeerde resultaten en verkeerde conclusies. Dit heeft grote gevolgen [1].
Er zijn vele parametrische en niet-parametrische modellen die toegepast kunnen worden. Een
parametrisch model is een verzameling van verdelingen die gebruik maakt van een eindig aantal
parameters. Een niet-parametrisch model maakt daarentegen gebruik van een oneindig aantal
parameters. Bij het prostaatkanker probleem dat aan het einde van dit verslag bekeken wordt
is geen vermoeden van een bepaald model. Daarom zal dat probleem geanalyseerd worden
door middel van niet-parametrische modellen. Om deze reden zal in het verslag alleen naar
niet-parametrische methoden gekeken worden.
In hoofdstuk 2 zal de survivalanalyse verder toegelicht worden. De belangrijke formules
worden gedefinieerd en de Kaplan-Meier schatter, de Nelson-Aalen schatter en de log-rank toets
en het Cox proportionele hazardmodel worden afgeleid. In hoofdstuk 3 wordt uitgelegd wat
competing risks precies zijn en worden twee competing risk modellen toegelicht. Deze modellen
zijn het oorzaak-specifieke hazardmodel en het Fine en Gray model.
In hoofdstuk 4 wordt de theorie toegepast op een prostaatkanker-dataset die beschikbaar
gesteld is door het ErasmusMC in Rotterdam. De dataset wordt geanalyseerd met de KaplanMeierschatter, het Fine en Gray model, de log-rank toets en het bekijken van het verschil in
percentages. In hoofdstuk 5 zijn twee datasets gesimuleerd. Bij dataset 1 heeft de ene groep
meer kans op de beoogde gebeurtenis dan de andere groep. Bij dataset 2 zijn de kansen op
gebeurtenissen gelijk. Beide datasets worden net als de prostaatkanker-dataset geanalyseerd
met de Kaplan-Meierschatter, het Fine en Gray model, de log-rank toets en het bekijken van
het verschil in percentages. In hoofdstuk 6 worden conclusies getrokken uit de resultaten die in
9
10
HOOFDSTUK 1. INLEIDING
hoofdstuk 4 en hoofdstuk 5 besproken zijn. Ten slotte worden in hoofdstuk 7 een aantal zaken
uit hoofdstuk 4 tot en met hoofdstuk 6 ter discussie gesteld.
Hoofdstuk 2
Survivalanalyse
Survivalanalyse is een methode om te analyseren hoe lang het duurt tot een gebeurtenis plaatsvindt. Survivalanalyse wordt veel toegepast in medische onderzoeken. Een paar voorbeelden
van onderzoeken waarbij survivalanalyse gebruikt wordt zijn
1. Een onderzoek naar de tijd tussen het ontstaan van een ziekte en het einde van de ziekte.
2. Een onderzoek naar de tijd tussen het genezen van een ziekte en een terugval, bijvoorbeeld
bij kanker.
2.1
Gecensureerde data
In survivalanalyse wordt vaak gewerkt met gecensureerde data. Data is gecensureerd als niet
bekend is wanneer een gebeurtenis precies plaatsgevonden heeft. Bij gecensureerde data wordt
onderscheid gemaakt tussen links-gecensureerde data, interval-gecensureerde data en rechtsgecensureerde data. Een gebeurtenis is links-gecensureerd als het tijdstip waarop de ziekte
begint onbekend is, doordat deze plaats gevonden heeft voordat het onderzoek begint [2]. Bij
interval-gecensureerde data is onbekend wanneer de gebeurtenis precies plaatsgevonden heeft,
het is echter wel bekend dat de gebeurtenis tussen twee bekende tijdstippen plaatsgevonden
heeft [2]. Een gebeurtenis is rechts-gecensureerd als de gebeurtenis nog niet plaatsgevonden
heeft op het laatst bekende tijdstip [2]. Daarnaast kan de data ook links afgeknot zijn. Een
gebeurtenis is links-afgeknot als de gebeurtenis al voor aanvang van het onderzoek is afgelopen,
waardoor de patient niet in het onderzoek meegenomen wordt [3]. In de rest van dit verslag zal
alleen ingegaan worden op rechts-gecensureerde data en links-afgeknotte data.
2.2
Functies
Bekend is dat de verdeling van een variabele X wordt gekarakteriseerd door een verdelingsfunctie.
De verdelingsfunctie wordt gekarakteriseerd door een aantal functies, die allemaal uit te drukken
zijn in de verdelingsfunctie. Eerst zullen een aantal variabelen gedefinieerd worden, daarna zullen
de functies gedefinieerd worden. Hierbij wordt gebruik gemaakt van [2] en [4].
Beschouw X1 , X2 , . . . , Xn uit een onbekende verdelingsfunctie F op [0, ∞), onafhankelijk
van C1 , C2 , . . . , Cn met bekende verdelingsfunctie H op [0, ∞). Hierbij is Xi het tijdstip waarop
de gebeurtenis die ondezocht wordt plaatsvindt en Ci het tijdstip van censuur voor persoon i.
Definieer vervolgens Ti = min{Xi , Ci }. De variabelen ∆i = 1{Xi ≤Ci } geeft aan of de gebeurtenis
plaatsvindt (∆i = 1) of dat de gebeurtenis rechts-gecensureerd is (∆i = 0).
11
12
HOOFDSTUK 2. SURVIVALANALYSE
De hazardfunctie λ(t) van een stochastische variabele T benadert de kans dat T in een kort
tijdsbestek na t valt, gegeven dat T ≥ t. De hazardfunctie wordt gegeven door
P(t ≤ T ≤ t + δt|T ≥ t)
δt
P(t ≤ T ≤ t + δt)
= lim
δt↓0
P(T ≥ t)δt
P(t ≤ T ≤ t + δt)
= lim
met F (t−) de linker limiet van F in t
δt↓0 (1 − F (t−))δt
1
S(t) − S(t + δt)
=
lim
S(t) δt↓0
δt
λ(t) = lim
δt↓0
=−
dS(t)
dt
S(t)
d log(S(t))
=−
dt
De discrete hazardfunctie λ(ti ) van een stochastische variabele T benadert de kans dat T op
tijdstip ti plaatsvindt, gegeven dat T > ti−1 . De discrete hazardfunctie wordt gegeven door
λ(ti ) = P(T ≤ ti |T > ti−1 )
(2.1)
De hazardfunctie wordt op een natuurlijke manier uitgedrukt in de survivalfunctie S(t),
welke wordt gegeven door
S(t) = P(T > t)
(2.2)
= 1 − F (t)
= e−Λ(t)
De survivalfunctie wordt op een natuurlijke manier uitgedrukt in de cumulatieve hazardfunctie Λ(t), welke gegeven wordt door
Z
t
Λ(t) =
λ(s)ds
(2.3)
0
Z
t
d log(S(s))
ds
ds
0
= − log(S(t))
=
−
Merk op dat als de tijd discreet is, dan volgt uit (2.1) dat de survivalfunctie als volgt
geschreven kan worden.
S(ti ) = P(T > ti )
= P(T > ti , T > ti−1 )
= P(T > ti |T > ti−1 )P(T > ti−1 )
= P(T > ti |T > ti−1 )S(ti−1 )
= (1 − P(T ≤ ti |T > ti−1 ))S(ti−1 )
= (1 − λ(ti ))S(ti−1 )
(2.4)
2.3. METHODEN
13
Vul vervolgens S(ti−1 ) = (1 − λ(ti−1 ))S(ti−2 ) in enzovoorts.
Hieruit volgt
S(ti ) =
i
Y
(1 − λ(tj ))
(2.5)
j=1
2.3
Methoden
In deze paragraaf wordt gekeken naar modellen voor survivalanalyse zonder covariaten en daarna
naar modellen voor survivalanalyse met covariaten. Covariaten zijn factoren waardoor (meetbare) verschillen tussen personen meegenomen kunnen worden in een model. Voorbeelden van
covariaten zijn leeftijd, geslacht, gewicht en leefgewoonten. Sommige covariaten kunnen interessant zijn in een onderzoek. Bijvoorbeeld als gekeken wordt naar het krijgen van longkanker, dan
kunnen leeftijd en roken interessante covariaten zijn.
Voor rechtsgecensureerde data wordt een onafhankelijkheidsaanname gedaan. Deze aanname
zegt dat voor ieder individu het tijdstip van censuur onafhankelijk is van het tijdstip waarop
de gebeurtenis plaatsvindt. Uit de onafhankelijkheidsaanname volgt dat de hazardfunctie voor
de individuen die gecensureerd zijn gelijk is aan de hazardfunctie voor de individuen waarbij de
gebeurtenis heeft plaatsgevonden [2].
2.3.1
Methoden zonder covariaten
Hieronder zullen de Kaplan-Meier schatter en de Nelson-Aalen schatter gedefinieerd worden [4].
Maar eerst zal onderstaande likelihood uitgewerkt worden [4].
Merk op dat de geobserveerde data gegeven kan worden door
Z1 , . . . , Zn met Zi = (Ti , ∆i )
= min(Xi , Ci ), 1{Xi ≤Ci }
Neem aan dat F en G dichtheden f en g hebben onder de Lebesgue-maat λ. De verdeling van
Z wordt gegeven door
P(T ≤ t, ∆ = 1) = P(T ≤ t, X ≤ C)
= P(X ≤ t, X ≤ C)
Z t Z ∞
=
f (x)h(c)dcdx
x=0 c=x
Z t
=
f (x)(1 − H(x))dx
x=0
Z t
=
f (x) − f (x)H(x)dx
x=0
Z t
= F (t) −
f (x)H(x)dx
x=0
Z t
= F (t) −
H(x)f (x)dx
x=0
14
HOOFDSTUK 2. SURVIVALANALYSE
en
P(T ≤ t, ∆ = 0) = P(T ≤ t, C ≤ X)
= P(C ≤ t, C ≤ X)
Z t Z ∞
f (x)h(c)dxdc
=
c=0 x=c
Z t
Z ∞
=
h(c)
f (x)dxdc
c=0
x=c
Z t
h(c)(1 − F (c))dc
=
c=0
Z t
h(c) − h(c)F (c)dc
=
c=0
Z t
= H(t) − F (t)H(t) +
H(x)f (x)dx
0
Z t
= H(t)(1 − F (t)) +
H(x)f (x)dx
0
Hieruit volgen de dichtheden van bovenstaande verdelingen
d
d
P(T ≤ t, ∆ = 1) =
dt
dt
!
Z
F (t) −
H(x)f (x)dx
[0,t]
= f (t) − H(t)f (t)
= f (t)(1 − H(t))
en
d
d
P(T ≤ t, ∆ = 0) =
dt
dt
!
Z
H(t)(1 − F (t)) +
H(x)f (x)dx
[0,t]
= h(t) − (h(t)F (t) + H(t)f (t)) + H(t)f (t)
= h(t) − h(t)F (t)
= H(t)(1 − F (t))
Hieruit volgt dat de dichtheid van Z kan worden gegeven door
g(t, ∆) = (f (t)(1 − H(t)))∆i (h(t)(1 − F (t)))1−∆i
In veel datasets komt het voor dat op tijdstip ti meerdere gebeurtenissen plaatsvinden.
Laat m ≤ n het aantal disjunctie tijdstippen ti in {t1 , t2 , . . . , tn } zijn. Laat u1 , u2 , . . . , um de
geordende disjuncte tijdstippen zijn.
Definieer nu
nj = #{i : ti = uj } het aantal keer dat uj in de geobserveerde tijdstippen zit
dj = #{i : ti = uj en di = 1} het aantal geobserveerde gebeurtenissen op uj
m
X
rj =
ni het totaal aantal geobserveerde individuen dat nog leeft op uj
i=j
(2.6)
(2.7)
2.3. METHODEN
15
Nu kan de log likelihood opgesteld worden en vervolgens worden omgeschreven met behulp
van (2.2) en (2.5)


m
Y
l(F ) = log  (f (tj )(1 − H(tj )))∆j (h(tj )(1 − F (tj )))1−∆j 
j=1
=
=
m
X
j=1
m
X
(∆j (log(f (tj )) + log(1 − H(tj ))) + (1 − ∆j )(log(h(tj )) + log(1 − F (tj ))))
(∆j log(f (tj )) + (1 − ∆j ) log(1 − F (tj )) + ∆j log(1 − H(tj ))
j=1
+ (1 − ∆j ) log(h(tj )))
m
(i) X
=
(∆j log(f (tj )) + (1 − ∆j ) log(1 − F (tj )) + c)
=
j=1
m
X
(dj log(λ(uj )) + (nj − dj ) log(S(uj )) + c)
j=1
=
=
m
X
j=1
m
X
(dj log(λ(uj )) − dj log(S(uj )) + nj log(S(uj )) + c)
(dj log(λ(uj )) − dj log(1 − λ(uj )) + nj log(S(uj )) + c)
j=1
=
m
X
!
j
Y
dj log(λ(uj )) − dj log(1 − λ(uj )) + nj log( (1 − λ(ui ))) + c
j=1
=
=
m
X
j=1
m
X
i=1
dj log(λ(uj )) − dj log(1 − λ(uj )) + nj
j
X
!
log(1 − λ(ui )) + c
i=1
(dj log(λ(uj )) − dj log(1 − λ(uj )) + c) + n1 log(1 − λ(u1 )) + n2 (log(1 − λ(u1 ))
j=1
+ log(1 − λ(u2 ))) + . . . + nm (log(1 − λ(u1 )) + log(1 − λ(u2 )) + . . .
+ log(1 − λ(um )))
m
X
=
(dj log(λ(uj )) − dj log(1 − λ(uj )) + c) + log(1 − λ(u1 ))(n1 + n2 + . . . + nm )
j=1
+ log(1 − λ(u2 ))(n2 + . . . + nm ) + . . . + log(1 − λ(um ))nm


m
m
X
X
dj log(λ(uj )) − dj log(1 − λ(uj )) + log(1 − λ(uj ))
=
ni + c
j=1
=
=
m
X
j=1
m
X
i=j
(dj log(λ(uj )) − dj log(1 − λ(uj )) + rj log(1 − λ(uj )) + c)
(dj log(λ(uj )) + (rj − dj ) log(1 − λ(uj )) + c)
j=1
Merk op dat (i) geldt doordat de log likelihood gemaximaliseerd wordt voor f . Daardoor zal
16
HOOFDSTUK 2. SURVIVALANALYSE
er niets veranderen aan de term ∆j log(1 − H(tj )) + (1 − ∆j ) log(h(tj )); deze term kan daarom
worden gezien als constante.
Deze log likelihood is maximaal voor de waarde
ˆ i ) = di
λ(u
ri
(2.8)
Kaplan-Meier schatter
De Kaplan-Meier schatter is een niet-parametrische maximum likelihood schatter voor de survivalfunctie bij rechtsgecensureerde data. De Kaplan-Meier schatter wordt ook wel de productlimiet schatter genoemd, en volgt uit het invullen van (2.8) in (2.5). Merk op dat de KaplanMeier schatter uitgaat van de onafhankelijkheidsaanname tussen X en C. Hieronder wordt de
Kaplan-Meier schatter gedefinieerd.
ˆ i) =
S(u
i
Y
ˆ i ))
(1 − λ(u
j=1
=
i Y
j=1
dj
1−
rj
Merk op dat (2.8) geen positieve waarde zal toewijzen aan tijdstippen waar alleen gecensureerde
waarnemingen zijn.
Toets voor de Kaplan-Meier schatter
Om er achter te komen of bijvoorbeeld een behandeling werkt, zijn er twee groepen nodig binnen
de data. Een onderzoeksgroep en een controlegroep. De onderzoeksgroep krijgt de behandeling
wel, de controlegroep niet. Voor deze twee groepen kan de Kaplan-Meier schatter bepaald
worden. Nu is de vraag hoe de verdelingen van de survivaltijd van de groepen zich tot elkaar
verhouden. Hiervoor worden de nulhypothese en de alternatieve hypothese opgesteld. Volgens
de nulhypothese zijn de verdelingen van de survivaltijden van beide groepen gelijk aan elkaar.
Volgens de alternatieve hypothese zijn de verdelingen van surivaltijd van beide groepen ongelijk
aan elkaar. Of volgens de eenzijdige alternatieve hypothese is de verdeling van de survivaltijd
van de onderzoeksgroep groter dan van de controlegroep. Verwerpen van de nulhypothese ten
gunste van de alternatieve hypothese gebeurt meestal als de p-waarde kleiner is dan 0.05 of 0.01.
De hypotheses kunnen getoetst worden door middel van de log-rank toets. De log-rank toets
vergelijkt het aantal geobserveerde gebeurtenissen op ui met het aantal verwachte gebeurtenissen
op ui . De log-rank toets kan gebruikt worden bij rechts-gecensureerde data. Hieronder wordt
de log-rank toets [4] gedefinieerd.
Laat u1 , u2 , . . . , um de geordende tijdstippen ti zijn met maximaal ´e´en gebeurtenis per tijdstip
en laat 1 ≤ i ≤ m.
Definieer nu
di0 het aantal geobserveerde gebeurtenissen in de controlegroep op tijdstip ui
di1 het aantal geobserveerde gebeurtenissen in de onderzoeksgroep op tijdstip ui
ri0 het aantal geobserveerde individuen dat nog in leven is in de controlegroep op tijdstip ui
ri1 het aantal geobserveerde individuen dat nog in leven is in de onderzoeksgroep op tijdstip ui
Merk op di = di0 + di1 met di gedefinieerd zoals in (2.6) en ri = ri0 + ri1 met ri gedefinieerd
zoals in (2.7).
2.3. METHODEN
17
De log-rank toetsingsgrootheid wordt nu gegeven door
Qw =
P
( m
w (d − e ))2
i=1
Pmi i1 2 i1
met wi > 0 gewichten,
i=1 wi vi
di
ei1 = ri1 ,
ri
ri di (ri − di )ri0
vi =
ri2 (ri − 1)
Nelson-Aalen schatter
De Nelson-Aalen schatter is een niet-parametrische schatter voor de cumulatieve hazardfunctie
bij rechtsgecensureerde data. De Nelson-Aalen schatter volgt uit het invullen van (2.8) in (2.3)
Merk op dat de Nelson-Aalen schatter uitgaat van de onafhankelijkheidsaanname tussen X en
C. Hieronder wordt de Nelson-Aalen schatter gedefinieerd.
Z t
ˆ
˜
Λ(t) =
λ(s)ds
0
X dj
=
rj
tj ≤t
Invullen van de Nelson-Aalen schatter in de formule voor de survivalfunctie (2.2) geeft tevens
een schatter voor de survivalfunctie.
˜
˜ = e−Λ(t)
S(t)
Y − dj
=
e rj
tj ≤t
2.3.2
Methode met covariaten
Als het effect van covariaten meegenomen moet worden in de analyse, zal gebruik gemaakt
moeten worden van een methode die rekening houdt met covariaten. Een veel gebruikte methode
hiervoor is het Cox proportionele hazardmodel. Het Cox proportionele hazardmodel is een
speciaal geval van het algemene proportionele hazardmodel. Laat c : (R) → [0, ∞) een functie
zijn met c(0) = 1. Het algemene proportionele hazardmodel wordt gegeven door
λ(t|x) = λ0 (t)c β T x met λ0 de baseline hazard,
β ∈ Rd de regressie parameter vector,
x ∈ Rd de covariaten vector
De baseline hazard is zo gedefinieerd dat als x = 0, dan is de hazardfunctie gelijk aan de baseline
hazard.
Als persoon x hazardfunctie h(t|x∗) heeft en persoon y hazardfunctie h(t|x), dan is de
hazardratio gedefinieerd door h(t|x∗)
h(t|x) . De hazardratio geeft de verhouding tussen de kans van
persoon x op de gebeurtenis op tijdstip t en de kans van persoon y op de gebeurtenis op tijdstip
t. Als h(t|x∗)
h(t|x) < 1, dan heeft persoon x een kleinere kans op de gebeurtenis op tijdstip t dan
h(t|x∗)
h(t|x) = 1, dan hebben persoon x en persoon
Als h(t|x∗)
h(t|x) > 1, dan heeft persoon x een grotere
persoon y. Als
y dezelfde kans op de gebeurtenis
op tijdstip t.
kans op de gebeurtenis op tijdstip
18
HOOFDSTUK 2. SURVIVALANALYSE
t dan persoon y. Het proportionele hazardmodel neemt aan dat de hazardratio constant is [5].
Oftewel de hazardfunctie van persoon x is evenredig met de hazardfunctie van persoon y. Dit
volgt wanneer de hazardratio uitgewerkt wordt.
h(t|x∗)
h0 (t)c(β T x∗)
=
h( t|x)
h0 (t)c(β T x)
=
c(β T x∗)
c(β T x)
Te zien is dat de aanname geldt zolang de covariaten onafhankelijk van de tijd zijn. Doordat de
hazardfuncties van verschillende personen evenredig zijn, zullen de hazards elkaar niet kruisen.
Een algemene regel is daarom ook: als de hazards kruisen, dan is het proportionele hazardmodel
niet toepasbaar. In de rest van het verslag zal deze aanname de proportionele hazardaanname
genoemd worden.
Voor de functie c wordt vaak gekozen voor c(y) = ey , dit is het Cox proportionele hazardmodel.
λ(t|x) = λ0 (t)eβ
Tx
met λ0 de baseline hazard,
(2.9)
d
β ∈ R de regressie parameter vector,
x ∈ Rd de covariaten vector
De baseline hazard is zo gedefinieerd dat als x = 0, dan is de hazardfunctie gelijk aan de baseline
hazard.
De survivalfunctie [6] wordt dan gegeven door
S(t|x) = eΛ0 (t)e
βT x
= S0e
βT x
met Λ0 de cumulatieve baseline hazard
met S0 de baseline survivalfunctie
De cumulatieve baseline hazard is zo gedefinieerd dat als x = 0, dan is de cumulatieve hazardfunctie gelijk aan de cumulatieve baseline hazard. De baseline survivalfunctie is zo gedefinieerd
dat als x = 0, dan is de survivalfunctie gelijk aan de baseline survivalfunctie.
Om λ(t|x) te schatten zullen eerst β en λ0 (t) geschat moeten worden. Om te beginnen
kan de regressie parameter vector β geschat worden door het maximaliseren van de parti¨ele log
likelihood voor β. De parti¨ele log likelihood voor β is het gedeelte van de log likelihood dat
afhangt van β. Hier volgt een schatter voor β uit. Deze schatter kan ingevuld worden in de log
likelihood. Nu kan door middel van het maximaliseren van de log likelihood ook λ0 (t) bepaald
worden. Hieronder zal eerst de parti¨ele log likelihood afgeleid worden, daarna zal ook de log
likelihood afgeleid worden [4].
2.3. METHODEN
19
Beschouw een steekproef X1 , X2 , . . . , Xk onafhankelijke continue random variabelen zo dat,
Xi de verdeling heeft met hazardfunctie λi . Laat δt > 0 heel klein zijn. Dan volgt onderstaande
conditionele kans.
pt,δt (i) = P(t ≤ Xi ≤ t + δt|Xj ≥ t∀j, t ≤ Xl ≤ t + δt voor precies ´e´en l)
=
P(t ≤ Xi ≤ t + δt, t ≤ Xl ≤ t + δt voor precies ´e´en l|Xj ≥ t∀j)
P(t ≤ Xl ≤ t + δt voor precies ´e´en l|Xj ≥ t∀j)
=
P(t ≤ Xi ≤ t + δt, Xl ≤ t + δt∀l 6= i|Xj ≥ t∀j)
P(t ≤ Xl ≤ t + δt voor precies ´e´en l|Xj ≥ t∀j)
=
P(t ≤ Xi ≤ t + δt, Xl ≤ t + δt∀l 6= i| ∩kl=1 Cl )
met Cl = {Xj ≥ t}
P(t ≤ Xl ≤ t + δt voor precies ´e´en l|Xj ≥ t∀j)
=
P(t ≤ Xi ≤ t + δt, Xl ≤ t + δt∀l 6= i, ∩kl=1 Cl )
met Cl = {Xj ≥ t}
P(∩kl=1 Cl )P(t ≤ Xl ≤ t + δt voor precies ´e´en l|Xj ≥ t∀j)
=
P(t ≤ Xi ≤ t + δt, ∩l6=i {Xl > t + δt}, ∩kl=1 Cl )
met Cl = {Xj ≥ t}
P(∩kl=1 Cl )P(t ≤ Xl ≤ t + δt voor precies ´e´en l|Xj ≥ t∀j)
=
{P(t ≤ Xi ≤ t + δt} ∩ Ci , ∩l6=i {Xl > t + δt} ∩ Cl )
P(kCi )P(∩l6=i Cl )P(t ≤ Xl ≤ t + δt voor precies ´e´en l|Xj ≥ t∀j)
met Cl = {Xj ≥ t}
=
=
P(t ≤ Xi ≤ t + δt|Ci )
Q
l6=i P(Xl
> t + δt|Cl )
met Cl = {Xj ≥ t}
P(t ≤ Xl ≤ t + δt voor precies ´e´en l|Xj ≥ t∀j)
Q
P(t ≤ Xi ≤ t + δt|Xi > t) l6=i P(Xl > t + δt|Xl > t)
P(t ≤ Xl ≤ t + δt voor precies ´e´en l|Xj ≥ t∀j)
λi (t)
= Pk
j=1 λj (t)
Nu wordt λi (t) vervangen door λi (t|x) en wordt aangenomen dat maar ´e´en gebeurtenis plaatsT
vindt per tijdstip. Laat Ri = {j : tj ≥ ti }. Invullen van λ(t|x) = λ0 (t)eβ x geeft dan
pt,δt (i) = P
λi (tl |xi )
j∈Ri λj (tl |xj )
T
λ0 (tl )eβ xi
=P
β T xj
j∈Ri λ0 (tl )e
=P
eβ
Tx
j∈Ri
i
eβ
Tx
(2.10)
j
20
HOOFDSTUK 2. SURVIVALANALYSE
Hieruit volgt de parti¨ele log likelihood.
l(β) = log
k
Y
eβ
P
i=1
=
k
X
∆i log
j∈Ri
=
!∆i
i
eβ
Tx
eβ
Tx
j
eβ
j∈Ri
!
i
P
i=1
k
X
Tx
Tx
j


∆i log eβ
Tx
i
− log 
i=1
=
k
X

X
eβ
Tx
j

j∈Ri


∆i β T xi − log 
i=1

X
eβ
Tx
j
(2.11)

j∈Ri
Hieruit volgt een schatter voor β. Deze schatter is consistent en asymptotisch normaal [7]. Merk
hierbij op dat het mogelijk is dat er meerdere gebeurtenissen plaatsvinden op tj . Als bekend
is in welke volgorde de gebeurtenissen plaatsgevonden hebben, dan kunnen de kansen van de
gebeurtenissen in die volgorde vermenigvuldigd worden in de likelihood van (2.10). In dit geval
zijn er geen problemen. Echter wanneer de volgorde van de gebeurtenissen onbekend is, zal de
uitkomst van de likelihood van (2.10) afhankelijk zijn van de volgorde die gekozen wordt. Dit
geeft wel een probleem, waarop hier niet verder in gegaan wordt [4].
Onder de proportionele hazardaanname kan de log likelihood van (λ0 , β) bepaald worden [5]
[4]. Aangezien er uit bovenstaande parti¨ele log likelihood een schatting voor β volgt, kan de
geschatte waarde voor β genomen worden. Onderstaande log likelihood hoeft daardoor alleen
nog te maximaliseren over λ0 .
n
Y
l(λ0 ) = log
!
λ(ti |xi ) S(ti |xi )
i=1
n Y
= log
∆i
λ0 (ti )eβxi
∆ i
!
eβxi
S0 (ti )
i=1
=
n X
βx
∆i log λ0 (ti )eβxi + log S0 (ti )e i
i=1
=
=
=
n X
i=1
n X
i=1
n X
∆i log(λ0 (ti )) + log eβxi
+ eβxi log(S0 (ti ))
∆i (log(λ0 (ti )) + βxi ) + eβxi log(S0 (ti ))
∆i log(λ0 (ti )) + ∆i βxi + eβxi log(S0 (ti ))
i=1
=
n
X
i=1

∆i log(λ0 (ti )) + ∆i βxi + eβxi log
i
Y

(1 − λj )
j=1
2.3. METHODEN
=
n
X
21

∆i log(λ0 (ti )) + ∆i βxi + eβxi
i=1
=
n
X
i
X

log(1 − λj )
j=1
(∆i log(λ0 (ti )) + ∆i βxi ) + eβ
Tx
1
log(1 − λ1 ) + eβ
Tx
2
(log(1 − λ1 )
i=1
T
+ log(1 − λ2 )) + . . . + eβ xn (log(1 − λ1 ) + log(1 − λ2 ) + . . . + log(1 − λn ))
n
X
T
T
T
=
(∆i log(λ0 (ti )) + ∆i βxi ) + log(1 − λ1 )(eβ x1 + eβ x2 + . . . + eβ xn )
i=1
T
T
+ log(1 − λ2 )(eβ x2 + . . . + eβ xn ) + . . . + log(1 − λn )eβ


n
n
X
X
T
∆i log(λ0 (ti )) + ∆i βxi + log(1 − λi )
=
e β xj 
i=1
Tx
n
j=i


n
n
X
X
T
(i)
∆i log(λ0 (ti )) + log(1 − λi )
=
eβ xj + c
i=1
j=i
Merk op dat (i) geldt doordat de log likelihood gemaximaliseerd wordt voor de term λ0 . Daardoor zal er niets veranderen aan ∆i βxi ; deze term kan daarom worden gezien als constante.
De log likelihood is maximaal voor de waarde
−1

ˆ 0 (ti ) = ∆i 1 +
λ
X
ˆT xj
eβ

j∈Ri
In plaats van de log likelihood uit te rekenen kan gebruik worden gemaakt van de Breslow
schatter. De Breslow schatter wordt gegeven door
−1

ˆ 0 (ti ) = ∆i 
λ
X
eβ
Tx
j

j∈Ri
Merk op dat de schatter die volgt uit de log likelihood en de Breslow schatter nagenoeg gelijk
zijn.
Nu β en λ0 bekend zijn kan λ(t|x) berekend worden met formule (2.9)
22
HOOFDSTUK 2. SURVIVALANALYSE
Hoofdstuk 3
Competing risk analyse
Zoals in hoofdstuk 2 vermeld is wordt in de survivalanalyse vaak gewerkt met rechts-gecensureerde data. Doodsoorzaken kunnen door verschillende redenen rechts gecensureerd zijn. De
redenen kunnen worden onderverdeeld in drie groepen.
1. Het onderzoek stopt op een vooraf bekend tijdstip
2. De persoon is uit het onderzoek gestapt
3. Er treedt een andere doodsoorzaak op waardoor de beoogde doodsoorzaak niet plaats kan
vinden
Een andere doodsoorzaak waardoor de beoogde doodsoorzaak niet plaats kan vinden heet een
competing risk. Een voorbeeld van een competing risk is sterfte door een hartaanval als sterfte
door kanker onderzocht wordt. Dit is een voorbeeld met twee mogelijke gebeurtenissen. Er zijn
echter vaak meer gebeurtenissen, die er voor kunnen zorgen dat een beoogde gebeurtenis niet
plaatsvindt. Al die gebeurtenissen zijn dan competing risks. Hieronder staat een schematische
weergave van een situatie met competing risks.
Doodsoorzaak 1
Beoogde doodsoorzaak
Levend
@
1 Doodsoorzaak 2
..
..
..
..
.
@
@
@
Competing risk
..
..
..
..
.
@
@
@
R Doodsoorzaak K
Competing risk
Figuur 3.1: Situatie met k mogelijke Doodsoorzaken
In de standaard survival modellen, zoals in hoofdstuk 2 beschreven, wordt de aanname
gemaakt dat het tijdstip van de beoogde doodsoorzaak onafhankelijk is van het tijdstip van
23
24
HOOFDSTUK 3. COMPETING RISK ANALYSE
censuur. Deze aanname is juist als het onderzoek op een vooraf gekozen tijdstip stopt. Echter in
het geval dat de persoon uit het onderzoek stapt of in het geval dat een competing risk optreedt
is deze aanname in het algemeen niet juist [2].
Als een persoon uit het onderzoek stapt kan dat zijn omdat de persoon zo ernstig ziek is
dat de persoon niet meer naar het ziekenhuis toe kan komen. In dit geval is de keuze om uit
het onderzoek te stappen dus afhankelijk van de kans op sterven aan de beoogde doodsoorzaak.
Merk op dat als het onderzoek stopt op een vooraf bekend tijdstip maar bij een persoon de
beoogde doodsoorzaak nog niet opgetreden is, dan kan deze persoon alsnog overlijden aan de
beoogde doodsoorzaak nadat het onderzoek afgelopen is. Maar als een persoon overlijdt aan een
competing risk, dan kan deze persoon niet meer overlijden aan de beoogde doodsoorzaak. Ook
gebeurt het vaak dat een behandeling voor een ziekte de kans op een competing risk verhoogt.
Hieruit volgt dat sterven aan een competing risk afhankelijk is van de kans op sterven aan de
beoogde doodsoorzaak.
Vanwege deze afhankelijkheid zullen bij data met competing risks andere modellen gebruikt
moeten worden dan de modellen die gebruikt worden bij standaard survival analyse. De KaplanMeier schatter en de Nelson-Aalen schatter, afgeleid in hoofdstuk 2, werken beide met de onafhankelijkheidsaanname tussen het tijdstip van de doodsoorzaak en het tijdstip van censuur,
en kunnen in situaties met competing risks dus niet zomaar gebruikt worden. Wanneer deze
modellen wel gebruikt worden, zal daar een overschating van de kans op sterven aan de beoogde
doodsoorzaak uit volgen. De overschatting zal groter zijn naar mate competing risks een grotere
rol spelen [2].
Voorbeelden van modellen die wel gebruikt kunnen worden in situaties met competing risks
zijn het oorzaak-specifieke hazard model en het model van Fine en Gray. Deze twee modellen
zullen in de volgende paragraaf besproken worden.
3.1
Methoden
Verschillende modellen kunnen gebruikt worden om competing risk analyse uit te voeren. De
modellen die hier bekeken worden zijn het oorzaak-specifieke hazard regressie model en het
model van Fine en Gray. Dit zijn de twee modellen die vaak gebruikt worden. Beide modellen
zijn niet-parametrische modellen. In de modellen is ∆ de doodsoorzaak die plaatsvindt en T
het tijdstip waarop doodsoorzaak ∆ plaats vindt.
3.1.1
Methode zonder covariaten
Het oorzaak-specifieke hazardmodel
Het oorzaak-specifieke hazardmodel maakt gebruik van de oorzaak-specifieke hazardfunctie. De
oorzaak-specifieke hazardfunctie geeft de kans dat een individu sterft door oorzaak k op het
tijdsinterval (t, t + δt), gegeven dat het individu nog leeft op tijdstip t. De oorzaak-specifieke
hazardfunctie voor doodsoorzaak k [2] wordt gegeven door
P(t ≤ T ≤ t + δt, ∆i = k|T ≥ t)
δt↓0
δt
λk (t) = lim
Hieronder staat een schematische weergave van de oorzaak-specifieke hazardfuncties.
3.1. METHODEN
25
Doodsoorzaak 1
λ1 (t)
Levend
@
Beoogde doodsoorzaak
1 Doodsoorzaak 2
λ2 (t)
..
..
..
..
.
@
@
@
Competing risk
..
..
..
..
.
@
λK (t)@
@
R Doodsoorzaak K
Competing risk
Figuur 3.2: oorzaak-specifieke hazardfuncties bij een competing risk situatie
De oorzaak-specifieke hazardfunctie bekijkt de kans op sterven door oorzaak k op tijdstip t.
Om te kans te berekenen op sterven door oorzaak k voor tijdstip t wordt de cumulatieve incidentiefunctie bepaald. De cumulatieve incidentiefunctie maakt gebruik van de oorzaak-specifieke
hazardfuncties. De cumulatieve incidentiefunctie van een stochastische variabele T benadert de
kans dat T voor tijdstip t valt en k de doodsoorzaak is. De cumulatieve incidentiefunctie voor
oorzaak k [2] wordt gegeven door
Ik (t) = P(T ≤ t, ∆i = k)
Z t
Z t
P
− K
Λl (t)
l=1
λk (s)e
λl (s)ds
=
ds met Λl (t) =
0
3.1.2
0
Methoden met covariaten
De functies die hierboven genoemd zijn houden geen rekening met verschil in leeftijd, geslacht,
gewicht, leefgewoonten of andere factoren waardoor (meetbare) verschillen tussen personen meegenomen kunnen worden in een model. Deze factoren worden covariaten genoemd. Sommige
covariaten kunnen interessant zijn voor een onderzoek. Bijvoorbeeld als gekeken wordt naar
het krijgen van longkanker, dan kunnen leeftijd en roken interessante covariaten zijn. De twee
meest gebruikte methoden zullen hieronder gedefinieerd en besproken worden. Dit zijn het
oorzaak-specifieke hazard regressie model en het Fine en Gray model.
26
HOOFDSTUK 3. COMPETING RISK ANALYSE
Het oorzaak-specifieke hazard regressie model
Het Cox proportionele hazard regressie model wordt net als in hoofdstuk 2 gebruikt om het
effect van covariaten te bekijken. Hieruit volgt de oorzaak-specifieke hazardfunctie [8].
P(t ≤ T ≤ t + δt, ∆i = k|T ≥ t, x)
δt↓0
δt
λk (t|x) = lim
T
= λk,0 (t)eβk x met λk,0 (t) de oorzaak-specifieke baseline hazard voor oorzaak k,
(3.1)
βk ∈ Rd de regressie parameter vector voor oorzaak k,
x ∈ Rd de covariaten vector
De oorzaak-specifieke baseline hazard is zo gedefinieerd dat als x = 0, dan is de oorzaak-specifieke
hazardfunctie gelijk aan de oorzaak-specifieke baseline hazard.
Om λk (t|x) te schatten met bovenstaande formule zullen eerst β en λk,0 (t) geschat moeten
worden. Om te beginnen kan de regressie parameter vector β geschat worden door het maximaliseren van de parti¨ele log likelihood. Hier volgt een schatter voor β uit. Deze schatter is
consistent en asymptotisch normaal [7]. Met behulp van de schatter voor β en de Breslow schatter kan ook λ0 (t) geschat worden. De parti¨ele log likelihood is exact gelijk aan de parti¨ele log
likelihood die afgeleid is in hoofdstuk 2.3.2, zie vergelijking (2.11). De Breslow schatter wordt
gegeven door
−1

ˆ k,0 (ti ) = ∆i 
λ
X
eβ
Tx
j

met Ri = {j : tj ≥ ti }
j∈Ri
De oorzaak-specifieke hazardfunctie bekijkt de kans op sterven aan oorzaak k op tijdstip t.
Om te kans te berekenen op sterven aan oorzaak k voor tijdstip t wordt de cumulatieve incidentiefunctie bepaalt. De cumulatieve incidentiefunctie maakt gebruik van de oorzaak-specifieke
hazardfuncties. De cumulatieve incidentiefunctie voor oorzaak k [8] wordt gegeven door
Ik (t|x) = P(T ≤ t, ∆i = k|x)
Z t
Z t
P
− K
Λk (t|x)
k=1
=
λk (s|x)e
ds met Λk (t|x) =
λk (s|x)ds
0
(3.2)
0
Merk op dat (3.2) niet alleen afhangt van de oorzaak-specifieke hazardfunctie van oorzaak
k, maar ook van de oorzaak-specifieke hazardfuncties van de andere oorzaken. Daarmee verliest
de cumulatieve incidentiefunctie de proportionele hazardaanname. Hieruit volgt dat er geen
directe relatie is tussen de oorzaak specifieke hazardfunctie en de cumulatieve incidentiefunctie.
Het gevolg is dat het effect van de covariaten op de cumulatieve incidentiefunctie niet bepaald
kan worden. Daardoor is de cumulatieve incidentiefunctie moeilijk te interpreteren en moet hier
voorzichtig mee omgegaan worden [2].
Het Fine en Gray model
Het Fine en Gray model is gemaakt om problemen van het oorzaak-specifieke hazardmodel
te voorkomen. Het model van Fine en Gray is een Cox proportioneel hazardmodel voor de
subdistributie, deze zal straks gedefinieerd worden. Het Fine en Gray model maakt gebruik van
de subdistributie hazardfunctie. De subdistributie hazardfunctie geeft de kans dat een individu
sterft door oorzaak k op het tijdsinterval (t, t + δt), gegeven dat het individu nog leeft op
3.1. METHODEN
27
tijdstip t of het individu overleden is aan een andere oorzaak dan oorzaak k. De subdistributie
hazardfunctie voor oorzaak k [7] wordt gegeven door
P(t ≤ T ≤ t + δt, ∆i = k|(T ≥ t) ∪ (T < t, ∆i 6= k)|x)
δt↓0
δt
λk (t|x) = lim
T
= λk,0 (t)eβk x met λk,0 (t)de subdistributie baseline hazard voor oorzaak k,
βk ∈ Rd de regressie parameter vector voor oorzaak k,
x ∈ Rd de covariaten vector
De subdistributie baseline hazard is zo gedefinieerd dat als x = 0, dan is de subdistributie
hazardfunctie gelijk aan de subdistributie baseline hazard.
Om λk (t|x) te schatten met bovenstaande formule zullen eerst β en λk,0 (t) geschat moeten
worden. Om te beginnen kan de regressie parameter vector β geschat worden door het maximaliseren van de parti¨ele log likelihood. Hier volgt een schatter voor β uit. Deze schatter is
consistent en asymptotisch normaal [7]. Met behulp van de schatter voor β en de Breslow schatter kan ook λ0 (t) geschat worden. De parti¨ele log likelihood is exact gelijk aan de parti¨ele log
likelihood die afgeleid is in hoofdstuk 2.3.2, zie vergelijking (2.11). De Breslow schatter wordt
gegeven door
−1

ˆ k,0 (ti ) = ∆i 
λ
X
eβ
Tx
j

met Ri = {j : tj ≥ ti ∪ (tj ≤ ti ∩ ∆i 6= k)}
j∈Ri
De oorzaak-specifieke hazardfunctie bekijkt de kans op sterven aan oorzaak k op tijdstip
t. Om de kans te berekenen op sterven aan oorzaak k voor tijdstip t wordt de cumulatieve
incidentiefunctie bepaald. De cumulatieve incidentiefunctie maakt gebruik van de oorzaakspecifieke hazardfuncties. De cumulatieve incidentiefunctie voor oorzaak k [7] wordt gegeven
door
Ik (t|x) = P(T ≤ t, ∆i = k|x)
= 1 − e−
Rt
0
λk (s|x)ds
Merk op dat de risicogroep voor oorzaak k behorend bij de subdistributie hazardfunctie onnatuurlijk is. Tot deze risicogroep behoren namelijk ook personen die voor tijdstip t al overleden
zijn aan een andere oorzaak dan oorzaak k. Daardoor is de subdistributie hazardfunctie moeilijk
te interpreteren en moet hier voorzichtig mee omgegaan worden.
Merk tevens op dat de cumulatieve incidentiefunctie alleen afhangt van de subdistributie
hazardfunctie van oorzaak k. Hieruit volgt dat er een directe relatie is tussen de subdistributie
hazardfunctie en de cumulatieve incidentiefunctie. Met als gevolg dat het effect van de covariaten
op de cumulatieve incidentiefunctie wel bepaald kan worden.
28
HOOFDSTUK 3. COMPETING RISK ANALYSE
Hoofdstuk 4
Analyse van prostaatkanker-data
4.1
Databeschrijving
De dataset die hier gebruikt wordt is verkregen voor een onderzoek naar de relevantie van
controle op prostaatkanker [9]. Voor dit onderzoek zijn 88283 mannen tussen 54 jaar en 75
jaar gevraagd deel te nemen. In totaal hebben 42376 mannen hier positief op gereageerd. Deze
mannen zijn aselect ingedeeld in twee groepen, de controlegroep en de onderzoeksgroep. Na
deze aselect verdeling behoren 21210 mannen tot de onderzoeksgroep en behoren 21166 mannen
tot de controlegroep.
Mannen in de onderzoeksgroep worden, indien zij nog geen prostaatkanker-diagnose hebben
en jonger dan 75 jaar zijn, iedere vier jaar getest op de waarde van het prostaat-specifiek antigeen
(PSA) in het bloed. Hoe hoger de PSA-waarde is, hoe groter de kans op prostaatkanker is. Indien
de PSA-waarde hoger is dan 3 ng/ml wordt een biopsie gedaan. Uit de biopsie blijkt of de man
prostaatkanker heeft. Als hieruit volgt dat de man prostaatkanker heeft wordt hij behandeld.
Als hieruit volgt dat de man geen prostaatkanker heeft wordt hij na een jaar opnieuw opgeroepen
voor de test. Mannen in de controlegroep worden niet getest, tenzij deze mannen met klachten
bij de huisarts komen. Mannen met prostaatkanker die overlijden worden onderzocht op de
doodsoorzaak, is de doodsoorzaak prostaatkanker of een competing risk.
Als een pati¨ent de diagnose prostaatkanker krijgt, wordt genoteerd wat zijn PSA waarde,
Ct-waarde, Cn-waarde, Cm-waarde en gleason-waarde zijn. De Ct-waarde geeft aan hoe ver
de tumor uitgebreid is. Deze waarde wordt aangegeven met T 1, T 2, T 3 en T 4; hoe hoger, hoe
verder de tumor uitgebreid is. De Cn-waarde geeft aan of de lymfeklieren aangetast zijn. N 0
betekent dat de lymfeklieren niet aangetast zijn, N 1 betekent dat de lymfeklieren wel aangetast
zijn. De Cm-waarde geeft aan of er uitzaaiingen zijn. M 0 betekent dat er geen uitzaaiingen zijn,
M 1 betekent dat er wel uitzaaiingen zijn. De Gleason-waarde geeft aan hoe agressief de tumor
is. Deze waarde ligt tussen de 6 en 10, hoe hoger deze waarde is, hoe agressiever de tumor is.
Het gemiddelde aantal volgjaren van de mannen is 12.8 jaar.
In de onderstaande analyse wordt gebruik gemaakt van 31752 mannen die aselect uit de
bovenstaande dataset beschikbaar gesteld zijn door het ErasmusMC in Rotterdam. In deze
overgebleven dataset zijn 14 mannen overleden, maar is onbekend wanneer zij precies overleden
zijn. Aangezien de datum van overlijden noodzakelijk is voor de analyse zijn deze mannen uit
de dataset verwijderd. Hierdoor zullen de resultaten die verkregen worden uit de analyse niet
geheel correct zijn. De dataset bestaat hierdoor nog uit 31738 mannen. Van deze mannen zitten
15871 mannen in de controlegroep en 15867 mannen in de onderzoeksgroep. De R-code waarin
de dataset bruikbaar gemaakt wordt voor de analyse is terug te vinden in de bijlage.
29
30
4.2
HOOFDSTUK 4. ANALYSE VAN PROSTAATKANKER-DATA
Resulaten
Zoals in de databeschrijving is vermeld, is de dataset verkregen om te onderzoeken of controleren
op prostaatkanker zin heeft. Om deze vraag in de conclusie te kunnen beantwoorden zal eerst
gekeken worden naar het verschil tussen de verhouding van het aantal overleden mannen per
groep. Vervolgens zal een analyse uitgevoerd worden door middel van de Kaplan-Meierschatter
en de log-rank toets. Daarna zal ook kort gekeken worden naar plots die gemaakt zijn met
behulp van het Fine en Gray model uit de competing risk analyse. Tot slot zal nog gekeken
worden naar het verschil tussen de verhoudingen van het aantal mannen dat overleden is aan
prostaatkanker per groep.
Bij alle toetsen die gebruikt worden, wordt uitgegaan van α = 0.05. Merk op dat als de
alternatieve hypothese tweezijdig is, dan moet worden uitgegaan van α = 0.1.
4.2.1
Aantal overleden mannen per groep
Omdat de controlegroep en de onderzoeksgroep aselect ingedeeld zijn wordt verwacht dat wel
of niet controleren op prostaatkanker het enige verschil tussen de groepen is. Een verschil
in de verhoudingen van het aantal overleden mannen tot en met 2010 per groep zou daarom
een indicatie kunnen geven of controleren op prostaatkanker zin heeft. Voor de controlegroep
is deze verhouding van het aantal overleden mannen tot en met 2010 ten opzichte van het
aantal mannen in de groep 0.30168. Voor de onderzoeksgroep is deze verhouding van het aantal
overleden mannen tot en met 2010 ten opzichte van het aantal mannen in de groep 0.30617.
Hieruit volgt dat in verhouding 0.00449 mannen minder overleden zijn in de controlegroep dan
in de onderzoeksgroep.
Om te bepalen of het verschil veroorzaakt wordt door het controleren op prostaatkanker,
worden de mannen in de controlegroep en de mannen in de onderzoeksgroep bij elkaar gevoegd en
vervolgens aselect in twee groepen verdeeld. In beide groepen zitten daardoor mannen die onder
controle geweest zijn en mannen die niet onder controle geweest zijn. Groep 1 bestaat, net als
de oorspronkelijke controlegroep, uit 15871 mannen. Groep 2 bestaat, net als de oorspronkelijke
onderzoeksgroep, uit 15867 mannen. Merk op dat de mannen opnieuw aselect verdeeld zijn in
de groepen. Daarom wordt er vanuit gegaan dat in beide groepen dezelfde verhouding mannen
overlijdt. Tevens kan verwacht worden dat in beide groepen hetzelfde aantal mannen zit van
gelijke leeftijd. Voor deze dataset met de opnieuw aselect ingedeelde groepen is het verschil in
de verhoudingen van het aantal overleden mannen tot en met 2010 ten opzichte van het aantal
mannen in de groep bepaald. Dit is 1000 keer herhaald. De R-code is terug te vinden in de
bijlage.
In onderstaande figuur zijn de verkregen resultaten door middel van een empirische verdelingsfunctie geplot. De rode punt die in de figuur te zien is, geeft het verschil in de verhoudingen
van de controlegroep en de onderzoeksgroep. Dit punt ligt op (−0.00449, 0.181). Hieruit volgt
dat 18.1 procent van de resultaten een lagere waarde heeft voor het verschil in de verhouding.
Merk op dat met behulp van deze plot de nulhypothese getoetst kan worden. De nulhypothese
zegt dat de verhouding van het aantal overleden mannen tot en met 2010 in de controlegroep
even groot is als in de onderzoeksgroep. De alternatieve hypothese zegt dat de verhouding in
de controlegroep groter is dan de verhouding in de onderzoeksgroep. Hieruit volgt dat er geen
significant verschil is tussen de verhoudingen van het aantal overleden mannen tot en met 2010
ten opzichte van het aantal mannen in de controlegroep en in de onderzoeksgroep.
4.2. RESULATEN
31
Figuur 4.1: empirische verdelingsfunctie op basis van de dataset voor het verschil in de verhouding van de overleden mannen tot en met 2010
4.2.2
Survivalanalyse
Als het zin heeft om te controleren op prostaatkanker, dan zal een duidelijk verschil zichtbaar
zijn tussen het percentage dat sterft aan prostaatkanker in de controlegroep en het percentage
dat sterft aan prostaatkanker in de onderzoeksgroep. Bekend is dat de mannen aselect zijn
ingedeeld in de controlegroep en in de onderzoeksgroep. Daarom wordt er vanuit gegaan dat in
beide groepen hetzelfde percentage sterft aan een competing risk. Tevens kan verwacht worden
dat in beide groepen hetzelfde aantal mannen zit van gelijke leeftijd. Hieruit volgt dat als het
zin heeft om te controleren op prostaatkanker, dan wordt verwacht dat een duidelijk verschil
zichtbaar zal zijn tussen het percentage dat sterft in de controlegroep en het percentage dat sterft
in de onderzoeksgroep. Omdat hierbij geen onderscheid gemaakt wordt in de doodsoorzaak,
kan hier gebruik gemaakt worden van een standaard survival methode. Er is gekozen voor de
Kaplan-Meier schatter. De R-code die gebruikt is om de plots te maken is terug te vinden in de
bijlage.
In de onderstaande figuur is de tijd in jaren uitgezet tegen de survivalkans. Daarbij worden
de gecensureerde tijdstippen aangegeven met kruisjes. Het onderzoek begint op t = 0. In de plot
is te zien dat de survivalfunctie van de controlegroep ongeveer gelijk is aan de survivalfunctie
van de onderzoeksgroep. Voor beide groepen is 17 jaar na aanvang van het onderzoek de
survivalkans ongeveer 0.6. Merk wel op dat de survivalfunctie van de controlegroep iets hoger
is dan de survivalfunctie van de onderzoeksgroep. Uitvoeren van de log-rank toets op deze
survivalfuncties geeft een toetsingsgrootheid van 0.438 met een p-waarde van 0.508. Merk op
dat de p-waarde bepaald is met behulp van de normale verdeling.
32
HOOFDSTUK 4. ANALYSE VAN PROSTAATKANKER-DATA
Figuur 4.2: Geschatte survivalfuncties van de controlegroep en van de onderzoeksgroep voor
prostaatkanker-dataset met behulp van de Kaplan-Meier schatter
In figuur 4.2 is voor alle mannen in het onderzoek de survivalfunctie te zien. In de groep
zitten mannen tussen de 54 jaar en de 75 jaar. Om een indruk te krijgen van het verschil
in de survivalfuncties per leeftijd is het volgende gedaan. Voor de controlegroep en voor de
onderzoeksgroep is per leeftijd bepaald na hoeveel jaar de survivalfunctie onder de 90 procent
komt.
In de onderstaande figuur zijn de verkregen leeftijden en de bijbehorende tijden tegen elkaar
uitgezet. Merk op dat in de figuur geen grote verschillen te zien zijn tussen de controlegroep en
de onderzoeksgroep. Wel is te zien dat de survivaltijd kleiner wordt naar mate de leeftijd hoger
is.
4.2. RESULATEN
33
Figuur 4.3: Geschatte tijdsduur per leeftijd tot een survivalkans van 90 procent voor
prostaatkanker-data, voor de controlegroep en voor de onderzoeksgroep met behulp van de
Kaplan-Meier schatter
Om te bepalen of het verschil in de survivalfuncties veroorzaakt wordt door het controleren op
prostaatkanker, worden de mannen uit de controlegroep en de mannen uit de onderzoeksgroep bij
elkaar gevoegd en vervolgens aselect in twee groepen verdeeld. In beide groepen zitten daardoor
mannen die onder controle geweest zijn en mannen die niet onder controle geweest zijn. Groep 1
bestaat, net als de oorspronkelijke controlegroep, uit 15871 mannen. Groep 2 bestaat, net als de
oorspronkelijke onderzoeksgroep, uit 15867 mannen. Merk op dat de mannen opnieuw aselect
verdeeld zijn in de groepen. Daarom wordt er vanuit gegaan dat in beide groepen hetzelfde
percentage mannen sterft. Tevens kan verwacht worden dat in beide groepen hetzelfde aantal
mannen zit van gelijke leeftijd. Voor deze dataset is opnieuw de Kaplan-Meier schatter bepaald.
Met het resultaat van de Kaplan-Meier schatter is vervolgens de log-rank toets weer uitgevoerd
en de toetsingsgrootheid opgeslagen. Dit is 1000 keer herhaald. De R-code is terug te vinden in
de bijlage.
In onderstaande figuur zijn de verkregen toetsingsgrootheden in een empirische verdelingsfunctie geplot. De rode punt die in de figuur te zien is, geeft de toetsingsgrootheid van de
log-tank toets die uitgevoerd is voor de controlegroep en de onderzoeksgroep. Dit punt ligt
op (0.438, 0.49). Hieruit volgt dat 49 procent van de toetsingsgrootheden een lagere waarde
heeft. Op basis van deze toets kan de nulhypothese, het heeft geen zin om te controleren op
prostaatkanker, niet verworpen worden.
34
HOOFDSTUK 4. ANALYSE VAN PROSTAATKANKER-DATA
Figuur 4.4: empirische verdelingsfunctie op basis van de dataset voor de log-rank toets
4.2.3
Competing risk analyse
Als het zin heeft om te controleren op prostaatkanker, dan zal een duidelijk verschil zichtbaar
zijn tussen het percentage mannen dat sterft door prostaatkanker in de controlegroep en het
percentage mannen dat sterft door prostaatkanker in de onderzoeksgroep. Het is tevens belangrijk om te bekijken wat er met de sterfte door competing risks gebeurt wanneer er gecontroleerd
wordt op prostaatkanker. Bekend is dat de mannen aselect zijn ingedeeld in de controlegroep
en in de onderzoeksgroep. Daarom wordt er vanuit gegaan dat in beide groepen hetzelfde aantal mannen zit van gelijke leeftijd. Omdat er bij deze vraag een onderscheid gemaakt wordt
in doodsoorzaak, moet hier gebruik gemaakt worden van een competing risk methode. Er is
gekozen voor het Fine en Gray model. De R-code die gebruikt is om de plots te maken is terug
te vinden in de bijlage.
In de onderstaande figuur zijn de cumulatieve incidentiefuncties gegeven voor de controlegroep en voor de onderzoeksgroep. In de figuur is te zien dat de cumulatieve incidentiefuncties
voor het overlijden aan prostaatkanker nagenoeg gelijk aan elkaar zijn. Ook de cumulatieve
incidentiefuncties voor de competing risks lijken erg op elkaar. Toch ligt de cumulatieve incidentiefunctie voor de competing risks van de onderzoeksgroep net iets hoger.
4.2. RESULATEN
35
Figuur 4.5: Cumulatieve incidentiefunctie voor prostaatkanker-dataset met behulp van de Fine
en Gray methode
Omdat er in figuur 4.5 geen verschil zichtbaar is tussen de cumulatieve incidentiefuncties
voor de beoogde gebeurtenis van de groepen, zijn deze hieronder uitvergroot weergegeven. Te
zien is dat op t = 6.5 het aantal sterfgevallen door prostaatkanker voor de beide groepen gelijk
is. Als t > 6.5 is het aantal sterfgevallen door prostaatkanker voor de controlegroep hoger dan
voor de onderzoeksgroep.
36
HOOFDSTUK 4. ANALYSE VAN PROSTAATKANKER-DATA
Figuur 4.6: Cumulatieve incidentiefunctie voor het overlijden aan prostaatkanker met behulp
van de Fine en Gray methode
Om figuur 4.6 te verduidelijken is tevens een tabel toegevoegd met de waarden van deze
twee cumulatieve incidentiefuncties op de tijdstippen: 5 jaar, 10 jaar en 15 jaar. Ook is het
verschil tussen de functies weergegeven, de waarde is positief als meer mannen gestorven zijn
aan prostaatkanker in de controle groep dan in de onderzoeksgroep en omgekeerd. Zoals te
zien is in de tabel, zijn de cumulatieve incidentiefuncties niet gelijk, maar zijn de verschillen
minimaal.
tijdstip in jaren
5
10
15
Controlegroep
0.000819104
0.004788608
0.009897190
Onderzoeksgroep
0.001386525
0.004033529
0.008075525
verschil
-0.000567
0.000755
0.00182
Tabel 4.1: Waarden van de cumulatieve incidentiefuncties voor de beoogde gebeurtenis op tijdstippen 5 jaar, 10 jaar en 15 jaar behorend bij de prostaatkanker-data
4.2.4
Aantal mannen dat overleden is aan prostaatkanker per groep
In figuur 4.5 is te zien dat de competing risks de doodsoorzaak domineren. Hierdoor is in
paragraaf 4.2.1 en in paragraaf 4.2.2 geen effect zichbaar bij het controleren op prostaatkanker.
In deze paragrafen wordt namelijk alleen gekeken naar alle overleden mannen per groep. Daarom
wordt de analyse van paragraaf 4.2.1 nu uitgevoerd voor het aantal mannen dat overleden is
aan prostaatkanker tot en met 2010 in verhouding tot het aantal mannen per groep. Voor de
controlegroep is deze verhouding 0.00895. Voor de onderzoeksgroep is deze verhouding 0.00725.
Hieruit volgt dat in verhouding 0.00170 mannen meer overleden zijn aan prostaatkanker in de
controlegroep dan in de onderzoeksgroep
Om te bepalen of het verschil veroorzaakt wordt door het controleren op prostaatkanker
4.2. RESULATEN
37
worden de mannen in de controlegroep en de mannen in de onderzoeksgroep bij elkaar gevoegd en
vervolgens aselect in twee groepen verdeeld. In beide groepen zitten daardoor mannen die onder
controle geweest zijn en mannen die niet onder controle geweest zijn. Groep 1 bestaat, net als
de oorspronkelijke controlegroep, uit 15871 mannen. Groep 2 bestaat, net als de oorspronkelijke
onderzoeksgroep, uit 15867 mannen. Merk op dat de mannen opnieuw aselect verdeeld zijn in
de groepen. Daarom wordt er vanuit gegaan dat in beide groepen dezelfde verhouding mannen
overlijdt. Tevens kan verwacht worden dat in beide groepen hetzelfde aantal mannen zit van
gelijke leeftijd. Voor deze dataset met opnieuw aselect ingedeelde groepen is het verschil in de
verhoudingen van het aantal mannen dat overleden is aan prostaatkanker tot en met 2010 ten
opzichte van het aantal mannen in de groep bepaald. Dit is 1000 keer herhaald. De R-code is
terug te vinden in de bijlage.
In onderstaande figuur zijn de verkregen resultaten door middel van een empirische verdelingsfunctie geplot. De rode punt die in de figuur te zien is, geeft het verschil in de verhoudingen
van de controlegroep en de onderzoeksgroep. Dit punt ligt op (0.00170, 0.968). Hieruit volgt
dat 3.2 procent van de resultaten een hogere waarde heeft voor het verschil in de verhoudingen.
Merk op dat met behulp van deze plot de nulhypothese getoetst kan worden. De nulhypothese
zegt dat de verhouding van het aantal mannen dat overleden is aan prostaakanker tot en met
2010 in de controlegroep even groot is als in de onderzoeksgroep. De alternatieve hypothese
zegt dat de verhouding in de controlegroep groter is dan de verhouding in de onderzoeksgroep.
Hieruit volgt dat er een significant verschil is tussen de verhoudingen van het aantal mannen
dat overleden is aan prostaatkanker tot en met 2010 ten opzichte van het aantal mannen in de
controlegroep en in de onderzoeksgroep. De nulhypothese wordt dus verworpen
Figuur 4.7: empirische verdelingsfunctie op basis van de dataset voor het verschil in de verhouding van de mannen die overleden zijn aan prostaatkanker tot en met 2010
38
HOOFDSTUK 4. ANALYSE VAN PROSTAATKANKER-DATA
Hoofdstuk 5
Simulatie
5.1
Simulatie beschrijving
Van de prostaatkanker-data is niet veel bekend. Zo is bijvoorbeeld onbekend welke hazardfuncties bij de data passen en welke covariaten in welke mate een rol spelen voor de survivalkans.
Dit kan hooguit geschat worden. Als gewerkt wordt met gesimuleerde data is dit wel bekend.
Voor het simuleren van de data is de R-code B. Haller en K. Ulm gebruikt [10] en aangepast.
De data wordt als volgt gesimuleerd. Om te beginnen worden functies gemaakt om twee oorzaakspecifieke hazardfuncties te bepalen die afhankelijk zijn van de tijd, regressieco¨effici¨enten en
covariaten. Hiervan is ´e´en hazardfunctie voor de beoogde gebeurtenis en ´e´en hazardfunctie voor
de competing risks. Ook worden de covariaten gesimuleerd. Daarna wordt per persoon 20 jaar
lang iedere maand bepaald of een gebeurtenis plaatsvindt, indien de persoon nog niet eerder
een gebeurtenis ondervonden heeft. Als een gebeurtenis plaatsvindt, wordt bepaald of dit de
beoogde gebeurtenis was of een competing risk. Tevens wordt het tijdstip van de gebeurtenis
genoteerd. Als na 20 jaar nog geen gebeurtenis plaatsgevonden heeft, is de gebeurtenis bij deze
persoon gecensureerd. Het tijdstip van censuur is dan dus 20 jaar. Nu is per persoon bepaald
op welk tijdstip de beoogde gebeurtenis, de competing risk of censuur plaatsvond. Vervolgens
wordt aan iedere persoon een willekeurige leeftijd toegekend. Deze leeftijd is minimaal 54 jaar en
maximaal 75 jaar. De gemiddelde leeftijd is 63 jaar en 7 maanden met een standaardafwijking
van 5 jaar en 7 maanden. Deze gegevens zijn gebasseerd op de prostaatkanker-data. Doordat de
leeftijd achteraf toegekend wordt, zijn de hazardfuncties niet afhankelijk van de leeftijd. Dit is
in de prostaatkanker-data wel het geval. De regressieco¨effici¨enten, covariaten en hazardfuncties
zijn eveneens niet hetzelfde als van de prostaatkanker-data. De R-code voor het simuleren van
de data is terug te vinden in de bijlage. Hier is tevens te vinden hoe de oorzaak-specifieke
hazardfuncties gedefinieerd zijn.
Met deze R-code zijn twee datasets gesimuleerd, elke dataset bestaat uit 30000 personen. Van
deze 30000 personen behoren 15000 personen tot de controlegroep en behoren 15000 personen tot
de onderzoeksgroep. In dataset 1 hebben de personen in de onderzoeksgroep slechts 50 procent
kans op de beoogde gebeurtenis in vergelijking met de personen in de controlegroep. In dataset
2 zijn de controlegroep en de onderzoeksgroep met dezelfde oorzaak-specifieke hazardfuncties
gesimuleerd.
5.2
Resultaten
Net als bij de resultaten van de data zal eerst gekeken worden naar het verschil tussen de
verhouding van het aantal overleden personen per groep. Vervolgens zal een analyse uitgevoerd
worden door middel van de Kaplan-Meierschatter en de log-rank toets. Daarna zal ook kort
39
40
HOOFDSTUK 5. SIMULATIE
gekeken worden naar plots die gemaakt zijn met behulp van het Fine en Gray model uit de
competing risk analyse. Tot slot zal nog gekeken worden naar het verschil tussen de verhoudingen
van het aantal personen dat overleden is aan de beoogde gebeurtenis per groep.
Dit wordt voor beide gesimuleerde datasets gedaan. Voor de dataset 1 geldt zoals eerder
gezegd; de personen in de onderzoeksgroep hebben slechts 50 procent kans op de beoogde gebeurtenis in vergelijking met de personen in de controlegroep. Omdat de groepen verder aselect
ingedeeld zijn, kan verwacht worden dat alleen het aantal personen dat overlijdt aan de beoogde
gebeurtenis anders is. Voor dataset 2 geldt zoals eerder gezegd; de personen in de controlegroep
en de onderzoeksgroep zijn met dezelfde oorzaak-specifieke hazardfuncties gesimuleerd. Ook
hier zijn de groepen aselect ingedeeld. Hierdoor kan verwacht worden dat het aantal personen
dat overlijdt in beide groepen gelijk is. Met deze kennis en achterliggende gedachten worden de
analyses in de rest van dit hoofdstuk gedaan.
Bij alle toetsen die gebruikt worden, wordt uitgegaan van α = 0.05. Merk op dat als de
alternatieve hypothese tweezijdig is, dan moet worden uitgegaan van α = 0.1.
5.2.1
Aantal overleden personen per groep
Dataset 1
Bij deze dataset kan verwacht worden dat er een verschil is in de verhoudingen van het aantal
overleden personen tot en met 2010 per groep. Voor de controlegroep is de verhouding van het
aantal overleden personen tot en met 2010 ten opzichte van het aantal personen in de groep
0.45580. Voor de onderzoeksgroep is deze verhouding van het aantal overleden personen tot en
met 2010 ten opzichte van het aantal personen in de groep 0.43467. Hieruit volgt dat in verhouding 0.02113 personen minder overleden zijn in de onderzoeksgroep dan in de controlegroep.
Om te bepalen of dit verschil veroorzaakt wordt door het verschil in de oorzaak-specifieke
hazardfunctie van de beoogde gebeurtenis, worden de personen in de controlegroep en de personen in de onderzoeksgroep bij elkaar gevoegd en vervolgens aselect in twee groepen verdeeld.
In beide groepen zitten daardoor personen uit de oorspronkelijke controlegroep en personen uit
de oorspronkelijke onderzoeksgroep. Groep 1 bestaat, net als de oorspronkelijke controlegroep,
uit 15000 personen. Groep 2 bestaat, net als de oorspronkelijke onderzoeksgroep, uit 15000
personen. Merk op dat de personen opnieuw aselect verdeeld zijn in de groepen. Daarom wordt
er vanuit gegaan dat in beide groepen dezelfde verhouding personen overlijdt. Voor deze dataset
met de opnieuw aselect ingedeelde groepen is het verschil in de verhoudingen van het aantal
overleden personen tot en met 2010 ten opzichte van het aantal personen in de groep bepaald.
Dit is 1000 keer herhaald. De R-code is terug te vinden in de bijlage.
In onderstaande figuur zijn de verkregen resultaten door middel van een empirische verdelingsfunctie geplot. De rode punt die in de figuur te zien is, geeft het verschil in de verhoudingen
van de controlegroep en de onderzoeksgroep. Dit punt ligt op (0.02113, 1). Hieruit volgt dat 0
procent van de resultaten een hogere waarde heeft voor het verschil in de verhouding. Merk op
dat met behulp van deze plot de nulhypothese getoetst kan worden. De nulhypothese zegt dat de
verhouding van het aantal overleden personen tot en met 2010 in de controlegroep even groot is
als in de onderzoeksgroep. De alternatieve hypothese zegt dat de verhouding in de controlegroep
groter is dan de verhouding in de onderzoeksgroep. Hieruit volgt dat er een significant verschil
is tussen de verhoudingen van het aantal overleden personen tot en met 2010 ten opzichte van
het aantal personen in de controlegroep en in de onderzoeksgroep.
5.2. RESULTATEN
41
Figuur 5.1: empirische verdelingsfunctie op basis van de gesimuleerde dataset 1 voor het verschil
in de verhouding van de overleden personen tot en met 2010
Dataset 2
Bij deze dataset kan verwacht worden dat de verhoudingen van het aantal overleden personen
tot en met 2010 per groep gelijk zijn. Voor de controlegroep is de verhouding van het aantal
overleden personen tot en met 2010 ten opzichte van het aantal personen in de groep 0.45340.
Voor de onderzoeksgroep is deze verhouding van het aantal overleden personen tot en met 2010
ten opzichte van het aantal personen in de groep 0.46120. Hieruit volgt dat in verhouding 0.0078
personen minder overleden zijn in de controlegroep dan in de onderzoeksgroep.
Om te bepalen of het verschil logisch is voor datasets waarin de oorzaak-specifieke hazardfuncties gelijk zijn, worden de personen in de controlegroep en de personen in de onderzoeksgroep bij elkaar gevoegd en vervolgens aselect in twee groepen verdeeld. In beide groepen zitten
daardoor personen uit de oorspronkelijke controlegroep en personen uit de oorspronkelijke onderzoeksgroep. Groep 1 bestaat, net als de oorspronkelijke controlegroep, uit 15000 personen.
Groep 2 bestaat, net als de oorspronkelijke onderzoeksgroep, uit 15000 personen. Merk op dat
de personen opnieuw aselect verdeeld zijn in de groepen. Daarom wordt er vanuit gegaan dat in
beide groepen dezelfde verhouding personen overlijdt. Voor deze dataset met de opnieuw aselect
ingedeelde groepen is het verschil in de verhoudingen van het aantal overleden personen tot en
met 2010 ten opzichte van het aantal personen in de groep bepaald. Dit is 1000 keer herhaald.
De R-code is terug te vinden in de bijlage.
In onderstaande figuur zijn de verkregen resultaten door middel van een empirische verdelingsfunctie geplot. De rode punt die in de figuur te zien is, geeft het verschil in de verhoudingen
van de controlegroep en de onderzoeksgroep. Dit punt ligt op (−0.0078, 0.102). Hieruit volgt
dat 1,02 procent van de resultaten een lagere waarde heeft voor het verschil in de verhouding.
Merk op dat met behulp van deze plot de nulhypothese getoetst kan worden. De nulhypothese
zegt dat de verhouding van het aantal overleden personen tot en met 2010 in de controlegroep
even groot is als in de onderzoeksgroep. De alternatieve hypothese zegt dat de verhouding in
42
HOOFDSTUK 5. SIMULATIE
de controlegroep groter is dan de verhouding in de onderzoeksgroep. Hieruit volgt dat er geen
significant verschil is tussen de verhoudingen van het aantal overleden personen tot en met 2010
ten opzichte van het aantal personen in de controlegroep en in de onderzoeksgroep.
Figuur 5.2: empirische verdelingsfunctie op basis van de gesimuleerde dataset 2 voor het verschil
in de verhouding van de overleden personen tot en met 2010
5.2.2
survivalanalyse
Dataset 1
In deze dataset kan verwacht worden dat er een duidelijk verschil zichtbaar is tussen het percentage dat sterft in de controlegroep en het percentage dat sterft in de onderzoeksgroep. Hierbij
wordt geen onderscheid gemaakt in de doodsoorzaak. Net als bij de prostaatkanker-data wordt
er daarom gebruik gemaakt van de Kaplan-Meier schatter. De R-code die gebruikt is om de
plots te maken is terug te vinden in de bijlage.
In de onderstaande figuur is de tijd in jaren uitgezet tegen de survivalkans. Daarbij worden de
gecensureerde tijdstippen aangegeven met kruisjes. In de plot is te zien dat de survivalfunctie van
de onderzoeksgroep hoger ligt dan de survivalfunctie van de controlegroep. Voor beide groepen
is in 20 jaar de surivalkans gedaald naar net geen 0.6. Uitvoeren van de log-rank toets op deze
surivalfuncties geeft een toetsingsgrootheid van 14.145 met p-waarde van 0.000169. Merk op dat
de p-waarde bepaald is met behulp van de normale verdeling.
5.2. RESULTATEN
43
Figuur 5.3: Geschatte survivalfuncties van de controlegroep en van de onderzoeksgroep voor
simulatie dataset 1 met behulp van de Kaplan-Meier schatter
In figuur 5.3 is voor alle personen in het onderzoek de survivalfunctie te zien. In de groep
zitten personen tussen de 54 jaar en de 75 jaar. Om een indruk te krijgen van het verschil
in de survivalfuncties per leeftijd is het volgende gedaan. Voor de controlegroep en voor de
onderzoeksgroep is per leeftijd bepaald na hoeveel jaar de surivalfunctie onder de 90 procent
komt.
In de onderstaande figuur zijn de verkregen leeftijden en de bijbehorende tijden tegen elkaar
uitgezet. Merk op dat de lijnen elkaar een aantal keer kruisen, maar de lijn van de onderzoeksgroep ligt gemiddeld iets hoger. Merk tevens op dat de tijd gelijk is voor de verschillende
leeftijden. Merk op dat de hazardfuncties, en dus ook de survivalfuncties, niet afhangen van de
leeftijd. Dit is terug te zien in de figuur doordat de lijnen geen helling hebben.
44
HOOFDSTUK 5. SIMULATIE
Figuur 5.4: Geschatte tijdsduur per leeftijd tot een survivalkans van 90 procent voor simulatie
dataset 1, voor de controlegroep en voor de onderzoeksgroep met behulp van de Kaplan-Meier
schatter
Om te bepalen of het verschil in de survivalfuncties veroorzaakt wordt door het verschil in
de oorzaak-specifieke hazardfuncties voor de beoogde gebeurtenis, worden de personen uit de
controlegroep en de personen uit de onderzoeksgroep bij elkaar gevoegd en vervolgens aselect
in twee groepen verdeeld. In beide groepen zitten daardoor personen uit de oorspronkelijke
controlegroep en personen uit de oorspronkelijke onderzoeksgroep. Groep 1 bestaat, net als
de controlegroep, uit 15000 personen. Groep 2 bestaat, net als de onderzoeksgroep, uit 15000
personen. Merk op dat de personen opnieuw aselect verdeeld zijn in de groepen. Daarom wordt
er vanuit gegaan dat in beide groepen hetzelfde percentage personen sterft. Voor deze dataset
is opnieuw de Kaplan-Meier schatter bepaald. Met het resultaat van de Kaplan-Meier schatter
is vervolgens de log-rank toets weer uitgevoerd en de toetsingsgrootheid opgeslagen. Dit is 1000
keer herhaald. De R-code is terug te vinden in de bijlage.
In onderstaande figuur zijn de verkregen toetsingsgrootheden in een empirische verdelingsfunctie geplot. De rode punt die in de figuur te zien is, geeft de toetsingsgrootheid van de
log-tank toets die uitgevoerd is voor de controlegroep en de onderzoeksgroep. Dit punt ligt op
(14.145, 1). Hieruit volgt dat 0 procent van de toetsingsgrootheden een hogere waarde heeft. Op
basis van deze toets kan de nulhypothese, het heeft geen zin om te controleren op prostaatkanker,
verworpen worden.
5.2. RESULTATEN
45
Figuur 5.5: empirische verdelingsfunctie op basis van gesimuleerde dataset 1 voor de log-rank
toets
Dataset 2
In deze dataset kan verwacht worden dat er geen verschil zichtbaar is tussen het percentage dat
sterft in de controlegroep en het percentage dat sterft in de onderzoeksgroep.
In de onderstaande figuur is de tijd in jaren uitgezet tegen de survivalkans. Daarbij worden de
gecensureerde tijdstippen aangegeven met kruisjes. In de plot is te zien dat de survivalfuncties
van de twee groepen ongeveer gelijk lopen. Voor beide groepen is in 20 jaar de surivalkans
gedaald naar net geen 0.6. Uitvoeren van de log-rank toets op deze surivalfuncties geeft een
toetsingsgrootheid van 1.573 met p-waarde van 0.21. Merk op dat de p-waarde bepaald is met
behulp van de normale verdeling.
46
HOOFDSTUK 5. SIMULATIE
Figuur 5.6: Geschatte survivalfuncties van de controlegroep en van de onderzoeksgroep voor
simulatie dataset 2 met behulp van de Kaplan-Meier schatter
In figuur 5.6 is voor alle personen in het onderzoek de survivalfunctie te zien. In de groep
zitten personen tussen de 54 jaar en de 75 jaar. Om een indruk te krijgen van het verschil
in de survivalfuncties per leeftijd is het volgende gedaan. Voor de controlegroep en voor de
onderzoeksgroep is per leeftijd bepaald na hoeveel jaar de surivalfunctie onder de 90 procent
komt.
In de onderstaande figuur zijn de verkregen leeftijden en de bijbehorende tijden tegen elkaar
uitgezet. Merk op dat in de figuur geen grote verschillen zichtbaar zijn tussen de controlegroep
en de onderzoeksgroep. Merk tevens op dat de tijd gelijk is voor de verschillende leeftijden.
Merk op dat de hazardfuncties, en dus ook de survivalfuncties, niet afhangen van de leeftijd.
Dit is terug te zien in de figuur doordat de lijnen geen helling hebben.
5.2. RESULTATEN
47
Figuur 5.7: Geschatte tijdsduur per leeftijd tot een survivalkans van 90 procent voor simulatie
dataset 2 voor de controlegroep en voor de onderzoeksgroep met behulp van de Kaplan-Meier
schatter
Om te bepalen of het verschil logisch is voor datasets waarin de oorzaak-specifieke hazardfuncties gelijk zijn, worden de personen uit de controlegroep en de personen uit de onderzoeksgroep bij elkaar gevoegd en vervolgens aselect in twee groepen verdeeld. In beide groepen zitten
daardoor personen uit de oorspronkelijke controlegroep en personen uit de oorspronkelijke onderzoeksgroep. Groep 1 bestaat, net als de controlegroep, uit 15000 personen. Groep 2 bestaat,
net als de onderzoeksgroep, uit 15000 personen. Merk op dat de personen opnieuw aselect
verdeeld zijn in de groepen. Daarom wordt er vanuit gegaan dat in beide groepen hetzelfde
percentage personen sterft. Voor deze dataset is opnieuw de Kaplan-Meier schatter bepaald.
Met het resultaat van de Kaplan-Meier schatter is vervolgens de log-rank toets weer uitgevoerd
en de toetsingsgrootheid opgeslagen. Dit is 1000 keer herhaald. De R-code is terug te vinden in
de bijlage.
In onderstaande figuur zijn de verkregen toetsingsgrootheden in een empirische verdelingsfunctie geplot. De rode punt die in de figuur te zien is, geeft de toetsingsgrootheid van de
log-tank toets die uitgevoerd is voor de controlegroep en de onderzoeksgroep. Dit punt ligt op
(1.5730, 0.791). Hieruit volgt dat 20.9 procent van de toetsingsgrootheden een hogere waarde
heeft. Op basis van deze toets kan de nulhypothese; het heeft geen zin om te controleren op
prostaatkanker, niet verworpen worden.
48
HOOFDSTUK 5. SIMULATIE
Figuur 5.8: empirische verdelingsfunctie op basis van dataset 2 voor de log-rank toets
5.2.3
Competing risk analyse
In tegenstelling tot paragraaf 5.2.1 en paragraaf 5.2.2, wordt nu wel onderscheid gemaakt in de
doodsoorzaak. Daarom moet een competing risk model gebruikt worden. Als competing risk
model is gekozen voor het Fine en Gray model. De R-code die gebruikt is om de plots te maken
is terug te vinden in de bijlage.
Dataset 1
In deze dataset kan verwacht worden dat er een duidelijk verschil zichtbaar is in de cumulatieve
incidentiefuncties voor prostaatkanker. Ook wordt verwacht dat er geen of weinig verschil is in
de cumulatieve incidentiefuncties van de competing risks.
In de onderstaande figuur zijn de cumulatieve incidentiefuncties gegeven voor de controlegroep en voor de onderzoeksgroep. In de figuur is te zien dat de cumulatieve incidentiefunctie
voor de beoogde gebeurtenis van de controlegroep een stuk hoger ligt dan de cumulatieve incidentiefunctie voor de beoogde gebeurtenis van de onderzoeksgroep. Ook de cumulatieve incidentiefuncties voor de competing risks verschillen van elkaar. In dit geval ligt de cumulatieve
incidentiefunctie van de controlegroep juist lager dan die van de onderzoeksgroep. Het verschil
tussen de cumulatieve incidentiefunctie voor de beoogde gebeurtenis is echter groter dan het
verschil tussen de cumulatieve incidentiefunctie voor competing risks.
5.2. RESULTATEN
49
Figuur 5.9: Cumulatieve incidentiefunctie voor gesimuleerde dataset 1 met behulp van de Fine
en Gray methode
Dataset 2
In deze dataset kan verwacht worden dat er geen of weinig verschil zichtbaar is in de cumulatieve
incidentiefuncties voor de beoogde gebeurtenis, en ook voor de competing risks.
In de onderstaande figuur zijn de cumulatieve incidentiefuncties gegeven voor de controlegroep en voor de onderzoeksgroep. In de figuur is te zien dat de cumulatieve incidentiefunctie
voor de beide gebeurtenissen ongeveer gelijk is tussen de groepen.
50
HOOFDSTUK 5. SIMULATIE
Figuur 5.10: Cumulatieve incidentiefunctie voor gesimuleerde dataset 2 met behulp van de Fine
en Gray methode
Omdat er in figuur 5.10 geen verschil zichtbaar is tussen de cumulatieve incidentiefuncties
voor de beoogde gebeurtenis van de groepen, zijn deze hieronder uitvergroot weergegeven. Te
zien is dat er inderdaad amper verschil zit tussen de cumulatieve incidentiefuncties.
Figuur 5.11: Cumulatieve incidentiefunctie voor het overlijden aan prostaatkanker met behulp
van de Fine en Gray methode
Hoofdstuk 6
Conclusie
Zoals in hoofdstuk 4 gezegd is, wordt er gezocht naar een antwoord op de vraag ’Heeft
zin om te controleren op prostaatkanker’. Om deze vraag te beantwoorden zullen eerst
aantal conclusies getrokken worden op basis van de prostaatkanker-data. Daarna zullen
aantal conclusies getrokken op basis van de simulatie. Het hoofdstuk wordt afgesloten met
slotconclusie waarin de eerdere conclusies gecombineerd worden.
6.1
het
een
een
een
Conclusies prostaatkanker-data
Als gekeken wordt naar het verschil in de verhoudingen van het aantal overleden mannen tot
en met 2010 ten opzichte van het aantal mannen in de groep, dan volgt uit figuur 4.1 een pwaarde van 0.181. Om deze p-waarde te krijgen is het verschil vergeleken met het verschil dat
volgt wanneer de groepen aselect ingedeeld worden. Op basis van dit resultaat volgt dat er
geen significant verschil is tussen de verhoudingen van het aantal overleden mannen tot en met
2010 ten opzichte van het aantal mannen in de controlegroep en in de onderzoeksgroep. Bij
deze resultaten kan de conclusie getrokken worden dat controleren op prostaatkanker geen effect
heeft.
In figuur 4.2 is te zien dat er weinig verschil is tussen de survivalfunctie van de controlegroep
en de survivalfunctie van onderzoeksgroep. De survivalfunctie van de controlegroep ligt echter
wel net iets hoger. Uit de log-rank toets volgt de toetsingsgrootheid 0.438 met p-waarde =
0.508. Op basis van deze resultaten kan de conclusie getrokken worden dat controleren op
prostaatkanker geen beter of slechter resultaat oplevert dan niet controleren op prostaatkanker.
Uit figuur 4.4 komt tevens een p-waarde van 0.49 als de toetsingsgrootheid vergeleken wordt met
de toetsingsgrootheden die volgen wanneer de groepen aselect ingedeeld worden. Dit bevestigt
bovenstaande conclusie.
In figuur 4.3 is te zien dat er weinig tot geen verschil is tussen de survivalkans van 0.9 op
iedere leeftijd voor de beide groepen. Hieruit kan de conclusie getrokken worden dat controleren
op prostaatkanker voor geen enkele leeftijd een beter of slechter resultaat oplevert.
Toch zijn de conclusies die hierboven getrokken worden niet wat verwacht wordt. In figuur
4.5, 4.6 en tabel 4.1 is te zien dat de kans op overlijden aan prostaatkanker in de controlegroep
niet gelijk is aan de kans op overlijden aan prostaatkanker in de onderzoeksgroep. De cumulatieve
incidentiefunctie voor de competing risks van de onderzoeksgroep is ongeveer gelijk aan de
cumulatieve incidentiefunctie voor de competing risks van de controlegroep.
Als nu gekeken wordt naar het verschil in de verhoudingen van het aantal mannen dat
overlijdt aan prostaatkanker tot en met 2010 ten opzichte van het aantal mannen in de groep,
dan volgt uit figuur 4.7 een p-waarde van 0.023. Om deze p-waarde te krijgen is het verschil
vergeleken met het verschil dat volgt wanneer de groepen aselect ingedeeld worden. Dit betekend
51
52
HOOFDSTUK 6. CONCLUSIE
dat 3.2 procent van de verschillen bij aselecte groepen groter is. Op basis van dit resultaat volgt
dat er een significant verschil is tussen de verhoudingen van het aantal mannen dat overlijdt
aan prostaatkanker tot en met 2010 ten opzichte van het aantal mannen per groep. Bij deze
resultaten kan de conclusie getrokken worden dat controleren op prostaatkanker toch wel effect
heeft.
6.2
Conclusies simulatie
Bij de simulatie zijn twee datasets gemaakt. In dataset 1 is wel een verschil tussen de controlegroep en de onderzoeksgroep. De onderzoeksgroep heeft namelijk maar 50 procent kans op de
beoogde gebeurtenis in vergelijking met de controlegroep. In dataset 2 is geen verschil tussen de
controlegroep en de onderzoeksgroep. Eerst zullen conclusies getrokken worden uit de resultaten
van dataset 1. Vervolgens zullen nog conclusies getrokken worden uit de resultaten van dataset
2. Bij de conclusies wordt tevens vermeld wat de conclusie zou zijn als de gesimuleerde datasets
de prostaatkanker-data zouden zijn.
6.2.1
Conclusies dataset 1
Als gekeken wordt naar het verschil in de verhoudingen van het aantal overleden personen tot
en met 2010 ten opzichte van het aantal personen in de groep, dan volgt uit figuur 5.1 een
p-waarde van 0. Om deze p-waarde te krijgen is het verschil vergeleken met het verschil dat
volgt wanneer de groepen aselect ingedeeld worden. Op basis van dit resultaat volgt dat er een
significant verschil is tussen de verhoudingen van het aantal overleden mannen tot en met 2010
ten opzichte van het aantal mannen in de controlegroep en in de onderzoeksgroep. Bij deze
resultaten kan de conclusie getrokken worden dat controleren op prostaatkanker effect heeft.
In figuur 5.3 is, zoals verwacht, een duidelijk verschil zichtbaar tussen de survivalfunctie
van de controlegroep en de survivalfunctie van de onderzoeksgroep. De survivalfunctie van de
onderzoeksgroep ligt hoger dan de survivalfunctie van de controlegroep. Uit de log-rank toets
volgt p-waarde = 0.000169, dit is ver onder de waarde α = 0.01. Op basis van deze resultaten
kan de conclusie getrokken worden dat er minder personen overlijden in de onderzoeksgroep
dan in de controlegroep. Met andere woorden: als deze dataset de prostaatkanker-data is, dan
heeft controleren op prostaatkanker zin. Uit figuur 5.5 komt tevens een p-waarde van 0 als de
toetsingsgrootheid vergeleken wordt met de toetsingsgrootheden die volgen wanneer de groepen
aselect ingedeeld worden. Dit bevestigt bovenstaande conclusie.
In figuur 5.9 is te zien dat de cumulatieve incidentiefuncties voor de beoogde gebeurtenis
verschillen. Dit is ook het geval voor de cumulatieve incidentiefuncties voor competing risks.
Voor de beoogde gebeurtenis ligt de cumulatieve incidentiefunctie van de beoogde gebeurtenis
lager. Voor competing risks ligt de cumulatieve incidentiefunctie van de beoogde gebeurtenis
juist hoger. Uit de survivalfunctie van figuur 5.3 en de cumulatieve incidentiefuncties van figuur
5.9 kan nu de conclusie getrokken worden dat de totale kans op overlijden in de onderzoeksgroep
kleiner is, dit komt doordat de kans op overlijden aan de beoogde gebeurtenis kleiner is. Maar
de kans op overlijden aan competing risks is wel groter.
In figuur 5.4 is te zien dat het langer duurt tot de survivalkans van 0.9 op iedere leeftijd voor
de onderzoeksgroep bereikt is dan voor de controlegroep. Ook dit resultaat komt overeen met
de conclusie dat er minder personen overlijden in de onderzoeksgroep dan in de controlegroep,
maar dat dit niet het geval is bij random groepen.
6.3. SLOTCONCLUSIE
6.2.2
53
Conclusies dataset 2
Als gekeken wordt naar het verschil in de verhoudingen van het aantal overleden personen tot
en met 2010 ten opzichte van het aantal personen in de groep. Dan volgt uit figuur 5.1 een
p-waarde van 0.102. Om deze p-waarde te krijgen is het verschil vergeleken met het verschil
dat volgt wanneer de groepen aselect ingedeeld worden. Op basis van dit resultaat volgt dat
er geen significant verschil is tussen de verhoudingen van het aantal overleden mannen tot en
met 2010 ten opzichte van het aantal mannen in de controlegroep en in de onderzoeksgroep. Bij
deze resultaten kan de conclusie getrokken worden dat controleren op prostaatkanker geen effect
heeft.
In figuur 5.6 is, zoals verwacht, bijna geen verschil te zien tussen de survivalfunctie van
de controlegroep en de survivalfunctie van de onderzoeksgroep. Uit de log-rank toets volgt pwaarde = 0.21. Op basis van deze resultaten kan de conclusie getrokken worden dat er evenveel
personen overlijden in de onderzoeksgroep als in de controlegroep. Met andere woorden: als deze
dataset de prostaatkanker-data is, dan heeft controleren op prostaatkanker geen invloed. Uit
figuur 5.8 komt tevens een p-waarde van 0.209 als de toetsingsgrootheid vergeleken wordt met
de toetsingsgrootheden die volgen wanneer de groepen aselect ingedeeld worden. Dit bevestigt
bovenstaande conclusie.
In figuur 5.10 is te zien dat de cumulatieve incidentiefuncties voor de beoogde gebeurtenis
gelijk zijn. Dit is ook het geval voor de cumulatieve incidentiefuncties voor competing risks.
Hieruit volgt dat de kans op sterven aan de beoogde gebeurtenis en de kans op sterven aan
competing risks tussen beide groepen niet verschillen. Dit bevestigt de conclusie dat er evenveel
personen overlijden in de onderzoeksgroep als in de controlegroep
In figuur 5.7 is te zien dat er weinig tot geen verschil is tussen de survivalkans van 0.9 op
iedere leeftijd voor de beide groepen. Ook dit resultaat komt overeen met de conclusie dat er
evenveel personen overlijden in de onderzoeksgroep als in de controlegroep, maar dat dit niet
het geval is bij random groepen.
6.3
Slotconclusie
Uit de resultaten en conclusies van de simulatie is te zien dat als er verschil is in het aantal
personen dat overlijdt aan de beoogde gebeurtenis, dit ook zichtbaar zal zijn in de plots en de
log-rank toets. Tevens is te zien dat als er geen verschil is in het aantal personen dat overlijdt
aan de beoogde gebeurtenis, en ook niet in het aantal personen dat overlijdt aan competing risks,
dit ook zichtbaar zal zijn in de plots en de log-rank toets. Bij de prostaatkanker-data volgt uit
globale analyses dat er geen effect zichtbaar is van het controleren op prostaatkanker. Zodra de
doodsoorzaken apart worden bekeken in de analyse blijkt dat het zin heeft om te controleren op
prostaatkanker.
54
HOOFDSTUK 6. CONCLUSIE
Hoofdstuk 7
Discussie
In hoofdstuk 4 is op prostaatkanker-data zowel een survivalanalyse als een competing risk analyse
toegepast. Daarbij was de vraag of controleren op prostaatkanker bij mannen tussen de leeftijd
van 54 jaar en 75 jaar zin heeft. In de conclusie zijn de resultaten die hieruit kwamen vergeleken
met de resultaten van de simulaties die in hoofdstuk 5 gedaan zijn. In hoofdstuk 6 is de conclusie
getrokken dat controleren op prostaatkanker zin heeft.
Hierbij moet opgemerkt worden dat de analyses uitgevoerd zijn op een random selectie van
31752 mannen uit een volledige dataset van 42376 mannen. Omdat het om een random selectie
gaat kan er vanuit gegaan worden dat dit geen invloed heeft op de resultaten. Dit kan echter nooit
met zekerheid gezegd worden. Daarbij zijn 14 mannen van de 31752 mannen niet meegenomen
in de analyses omdat van deze mannen geen datum van overlijden bekend is, terwijl zij wel
overleden zijn. Ook dit kan van invloed zijn op de resultaten.
Bij de analyse is naar survivalanalyse gekeken en naar een enkele plot vanuit competing risks
analyse. De conclusies zouden nauwkeuriger getrokken kunnen worden als er verder gekeken
wordt naar de competing risk analyse. Een toets bij de cumulatieve incidentiefuncties zou al
meer duidelijk kunnen maken. Het is nu bijvoorbeeld onbekend of het onderlinge verschil dat te
zien is tussen de cumulatieve incidentiefuncties voor de competing risks significant is. Het zou
namelijk een opmerkelijke uitkomst zijn als er gezegd kan worden dat er een significant verschil
is in de cumulatieve incidentiefuncties voor de competing risks, maar dat er geen significant
verschil is in de cumulatieve incidentiefuncties voor het overlijden aan prostaatkanker. In dat
geval zou het betekenen dat controleren op prostaatkanker een negatief effect heeft.
Binnen de analyse zijn leeftijd, PSA-waarde en andere covariaten niet meegenomen. Als
deze covariaten wel meegenomen worden in de analyse is het mogelijk dat er andere conclusies
getrokken moeten worden.
Verder zijn er nog een aantal interessante vragen die nog beantwoord kunnen worden door
middel van literatuurstudie. Ten eerste kan er nog gekeken worden naar andere competing risk
methoden. Voorbeelden hiervan zijn mixture models, vertical modelling en pseudo-observaties.
Ook is het interessant om uit te zoeken wat voor toetsen gebruikt kunnen worden bij competing
risk modellen en hoe deze werken.
55
56
HOOFDSTUK 7. DISCUSSIE
Bijlage A
Functies
# Functie voor het maken van een dataframe met daarin wanneer een bepaalde
# survivalkans gepasseerd wordt per leeftijd.
#
# Parameter dataset: Dataset waar de kolommen ftime en Death2010 uit komen
# Parameter groep: 0 (= controle groep) of 1 (= interventiegroep)
# Parameter kans: de surivalkans waarnaar gekeken wordt
tijdenSurvivalKans <- function(dataset, groep, kans){
# Ondergrens leeftijd
ondergrens <- 54
# Bovengrens leeftijd
bovengrens <- 75
# N: het aantal leeftijden
N <- bovengrens - ondergrens + 1
# Leeg dataframe aanmaken met: N rijen, kolom leeftijd, kolom tijd
dataTijden <- data.frame(Leeftijd = numeric(N), Tijd = numeric(N))
for(i in ondergrens:bovengrens){
# Kaplan-Meier schatter bepalen per leeftijd
surv <- survfit(Surv(ftime, Death2010)~group,
data = subset(dataset,
Age_at_rand >= i & Age_at_rand < (i+1)),
conf.type="none")
# Aanroepen functie tijdSurvivalKans
tijd <- tijdSurvivalKans(groep, surv, kans)
# Rij vullen in dataframe
dataTijden[i-ondergrens+1,] <- c(i, tijd)
}
return(dataTijden)
}
57
58
BIJLAGE A. FUNCTIES
# Functie om te bepalen op welk tijdstip de survivalkans onder
# een bepaalde waarde gaat
#
# Parameter groep: 0 (= controle groep) of 1 (= interventiegroep)
# Parameter survival: een survfit object
# Parameter kans: de surivalkans waarnaar gekeken wordt
tijdSurvivalKans <- function(groep, survival, kans){
index <- 1
rows <- NROW(survival[groep]$time)
# Bepaalt wanneer de survival kans lager is dan kans
while(survival[groep]$surv[index] > kans){
if(index > rows){
return(NA)
}
index <- index + 1
}
return(survival[groep]$time[index])
}
# Functie om de dataset te simuleren per groep
#
# Parameter aantalPersonen: aantal personen dat gesimuleerd wordt
# Parameter groep: 0 (= controle groep) of 1 (= interventiegroep)
# Parameter percentage: percentage van de oorspronkelijke kans op de
#
beoogde gebeurtenis
simuleren <- function(aantalPersonen, groep, percentage){
n <- aantalPersonen
# Aanmaken van gamma functie voor in lambda2_0
gamma_0 <- function(t) 0.001 * exp(-0.001*t/log(1.3))
# Definieren van oorzaak-specifieke baseline hazard voor de beoogde
# gebeurtenis
lambda1_0 <- function(t) percentage*0.001
# Definieren van de oorzaak-specifieke baseline hazard voor de competing risk
lambda2_0 <- function(t) gamma_0(t) - lambda1_0(t) grad(gamma_0,t)/gamma_0(t) +
grad(lambda1_0,t)/lambda1_0(t)
# Definieren van de regressie coefficienten beta_k (=log-hazard ratios)
xi <- c(log(1.15),log(0.9),log(2))
# Aanmaken van gamma functie voor in lambda2_x
gamma_x <- function(t,xi,X) gamma_0(t)*exp(sum(xi*X))
59
# Definieren van de oorzaak-specifieke hazardfunctie voor de beoogde
# gebeurtenis
lambda1_x <- function(t,X) 0.5*(lambda1_0(t)*exp(xi[1]*X[1]*exp(-0.0005*t)
+ xi[2]*X[2]*exp(-0.0005*t)
+ xi[3]*X[3]*exp(-0.0005*t)))
# Definieren van de oorzaak-specifieke hazardfunctie voor de competing risk
lambda2_x <- function(t,xi,X) 0.5*(gamma_x(t,xi,X) - lambda1_x(t,X)
- grad(gamma_x,t,xi=xi,X=X)/gamma_x(t,xi,X)
+ grad(lambda1_x,t,X=X)/lambda1_x(t,X))
# Genereren van covariaten
# X2 heeft een aselecte waarde tussen 0 en 10
X1 <- runif(n,0,3)
X2 <- pmax(pmin(rnorm(n,5,1),10),0)
X3 <- sample(0:1,n,repl=T)
XX <- cbind(X1,X2,X3)
# Vectoren om resultaten op te slaan vanuit onderstaande for-loop
Evstat <- rep(0,n)
Evtime <- c(0,n)
Event <- c()
# Bepalen van gebeurtenis, en het moment waarop die plaatsvindt per persoon
for(j in 1:n)
{
Ti <- 0
Event[j] <- 0
# while-loop loopt zolang er nog geen 20 jaar voorbij zijn of nog geen
# gebeurtenis plaatsgevonden heeft.
while(Ti< 20*12 && Event[j] == 0)
{
Ti <- Ti+1
# Kans op een gebeurtenis voor persoon j op tijdstip Ti
Prob <- lambda1_x(t=Ti,X=XX[j,]) + lambda2_x(t=Ti,xi=xi,X=XX[j,])
# Bepaalt of een gebeurtenis plaatsvindt
Event[j] <- sample(0:1,1,prob=c(1-Prob,Prob))
# Als een gebeurtenis plaatsvindt, bepaal welke gebeurtenis en tijdstip
if(Event[j]==1){
Evstat[j] <- sample(1:2,1,prob=c(
lambda1_x(t=Ti,X=XX[j,]),
lambda2_x(t=Ti,xi=xi,X=XX[j,])))
Evtime[j] <- Ti/12
}
else
Evtime[j] <- Ti/12
}
}
60
BIJLAGE A. FUNCTIES
# Kent aselect leeftijden tussen de 54 en 76 toe
# Gemiddelde leeftijd is 63.61, sd is 5.61
leeftijd <- pmax(pmin(rnorm(n,63.61,5.607905),75.99),54)
# Maak van de data een dataframe
simulatie <- data.frame(group = groep, Age_at_rand = leeftijd,
Death2010 = Event, gebeurtenis= Evstat,
ftime = Evtime)
# Geeft de dataframe terug
return(simulatie)
}
# Functie die een aantal toetsen uitvoert op de meegegeven data.
# De log-rank toets wordt uitgevoerd en het verschil in verhoudingen
# wordt uitgerekend.
#
# Parameter data: De data waar de toetsen op worden uitgevoerd.
# Parameter n0: het aantal personen in de controlegroep
# Parametet n1: het aantal personen in de onderzoeksgroep
toetsen <- function(data, n0, n1) {
# Vectoren om resultaten op te slaan vanuit onderstaande for-loop
logRank <- c()
verschilVerhoudingDeath <- c()
verschilVerhoudingPCAdeath <- c()
groepData <- c(rep(0, n0),rep(1, n1))
for(i in 1:10) {
# toets data maken waarin randGroup zit
tData <- data
tData$randGroup <- sample(groepData, (n0 + n1), replace = FALSE)
# Log-rank toets
logRankToets <- survdiff(Surv(ftime, Death2010)~randGroup,
data = tData, rho = 0)
logRank[i] <- logRankToets$chisq
# Opsplitsen van de controlegroep en de onderzoeksgroep
g1 <- subset(tData, randGroup == 0)
g2 <- subset(tData, randGroup == 1)
# Dataset met aantal personen dat overleden is per groep
death_g1 <- subset(g1, Death2010 == 1)
death_g2 <- subset(g2, Death2010 == 1)
61
# Dataset met aantal personen dat overleden is aan prostaatkanker per groep
PCA_death_g1 <- subset(g1, status == 1)
PCA_death_g2 <- subset(g2, status == 1)
# Aantal rijen in groep
ng1 <- nrow(g1)
ng2 <- nrow(g2)
# Bepalen van het verschil in de verhoudingen
verschilVerhoudingDeath[i] <- (nrow(death_g1) / ng1) (nrow(death_g2) / ng2)
verschilVerhoudingPCAdeath[i] <- (nrow(PCA_death_g1) / ng1) (nrow(PCA_death_g2) / ng2)
}
toetsRes <- data.frame(logRank, verschilVerhoudingDeath,
verschilVerhoudingPCAdeath)
return(toetsRes)
}
# Functie die de resultaten uit de functie toetsen plot dmv een empirische
# verdelingsfunctie.
#
# Parameter data: De data waar de toetsen op worden uitgevoerd.
# Parameter puntX: De toetsingsgrootheid die volgt uit de oorspronkelijke
#
prostaatkanker-dataset.
plotToetsData <-function(data, puntX) {
# Maakt een functie voor de empirische verdelingsfunctie
f <- ecdf(data)
# Plot de empirische verdelingsfunctie
plot(f, verticals = TRUE, do.points = FALSE, main = "")
# Plot de toetsingsgrootheid van de oorspronkelijke prostaatkanker-dataset
points(puntX, f(puntX), pch = 19, col = "red")
ablineclip(h = f(puntX), v = puntX, lty = 2, col = "red",
x2 = puntX, y2 = f(puntX))
}
62
BIJLAGE A. FUNCTIES
Bijlage B
Dataset
#### Dataset inlezen ####
data <- data.frame(as.data.set(
spss.system.file(’Data/Dataset/Competing_risk.sav’)))
#### Data 1 ####
# Dataset bruikbaar maken.
data1 <- data
# Dataset waarin de datum goed staat (niet numeriek)
data1$rand_date <- as.Date(ISOdate(1582,10,14) + data1$rand_date)
data1$deathdate <- as.Date(ISOdate(1582,10,14) + data1$deathdate)
data1$diag_date <- as.Date(ISOdate(1582,10,14) + data1$diag_date)
# Onbekende gleason waarde wordt aangegeven met NA
data1$ct <- replace(data1$ct, data1$ct == 999 | data1$ct == "", NA)
data1$cn <- replace(data1$cn, data1$cn == 99 | data1$cn == "", NA)
data1$cm <- replace(data1$cm, data1$cm == 99 | data1$cm == "", NA)
data1$gleason <- replace(data1$gleason, data1$gleason == 99 |
data1$gleason == "", NA)
# Deathdate en diag_date na 2010 wordt op NA gezet
data1$deathdate <- replace(data1$deathdate,
data1$deathdate > as.Date("2010-12-31", "%Y-%m-%d"),
NA)
data1$diag_date <- replace(data1$diag_date,
data1$diag_date > as.Date("2010-12-31", "%Y-%m-%d"),
NA)
# Toevoegen van kolom ftime aan data1
# ftime is het aantal dagen dat een patient gevolgd wordt vanaf het moment van
# de diagnose
# Hierbij is aangenomen dat patienten die nog leven gevolgt zijn tot 2010-12-31
ftime <- c()
cdate <- as.Date("2010-12-31", "%Y-%m-%d")
for(i in 1:nrow(data1)){
if(is.na(data1$deathdate[i])){
63
64
BIJLAGE B. DATASET
ftime[i] = cdate - data1$rand_date[i]
}
else {
ftime[i] = data1$deathdate[i] - data1$rand_date[i]
}
}
data1$ftime <- unlist(ftime/365.25)
# Aanmaken leeftijd waarop de personen uit het onderzoek verlaten door
# overlijden
# of andere redenen
Age_at_end <- c()
for(i in 1:nrow(data1)){
Age_at_end[i] = data1$Age_at_rand[i] + data1$ftime[i]
}
data1$Age_at_end <- unlist(Age_at_end)
# Toevoegen van kolom status van data1
# Status is 0 als patient nog leeft in 2010 (rechts gecensureerd).
# Status is 1 als patient aan prostaat kanker overleden is in 2010.
# Status is 2 als patient aan competing risk overleden is in 2010.
# Hierbij is aangenomen dat patienten die na 2010 overleden zijn
#status 0 hebben.
status <- c()
for(i in 1:nrow(data1)){
if(data1$Death2010[i] == 0){
status[i] = 0
}
else if(data1$Death2010[i] == 1 && data1$PCAdeath2010[i] == 1){
status[i] = 1
}
else{
status[i] = 2
}
}
data1$status <- unlist(status)
#### Data 2 ####
# Data waarbij personen wel overleden zijn, maar datum van overlijden onbekend
# is verwijderen
data2 <- data1
data2 <- subset(data2, !(is.na(deathdate) == TRUE & Death2010 == 1))
save(data2, file = "Data/Dataset/data2.R")
Bijlage C
Data simuleren
source(file="functies/bep_functies.r")
#### Simuleren data waarbij controle geen zin heeft ####
aantalPersonen <- 15000
# Simuleer de controle groep
controle <- simuleren(aantalPersonen, 0, 1)
# Simuleer de interventiegroep
interventie <- simuleren(aantalPersonen, 1, 1)
# Voeg de controlegroep en interventiegroep samen in 1 dataframe
dataSimulatie <- rbind(controle, interventie)
# De dataframe opslaan
save(dataSimulatie, file = "Simulatie/dataSimulatie.Rdata")
#### Simuleren data waarbij controle zin heeft ####
aantalPersonen <- 15000
# Simuleer de controle groep
controle <- simuleren(aantalPersonen, 0, 1)
# Simuleer de interventiegroep
# Deze groep heeft een lagere kans op sterven aan beoogde gebeurtenis
interventie <- simuleren(aantalPersonen, 1, 0.5)
# Voeg de controlegroep en interventiegroep samen in 1 dataframe
dataSimulatieControle <- rbind(controle, interventie)
# De dataframe opslaan
save(dataSimulatieControle, file = "Simulatie/dataSimulatieControle.Rdata")
65
66
BIJLAGE C. DATA SIMULEREN
Bijlage D
Analyse met verhoudingen
D.1
Prostaatkanker-dataset
source(file="functies/bep_functies.r")
load("Data/Dataset/data2.R")
#### Bepalen verhoudingen verschil ####
controle <- subset(data2, group == "control")
onderzoeks <- subset(data2, group == "intervention")
#Dataset met aantal personen dat overleden is per groep
Death2010Controle <- subset(controle, Death2010 == 1)
Death2010Onderzoeks <- subset(onderzoeks, Death2010 == 1)
#Dataset met aantal personen dat overleden is aan prostaatkanker per groep
PCAdeath2010Controle <- subset(controle, PCAdeath2010 == 1)
PCAdeath2010Onderzoeks <- subset(onderzoeks, PCAdeath2010 == 1)
# Bepaal verschil in de verhoudingen
verschilVerhoudingDeathData <- nrow(Death2010Controle)/nrow(controle) nrow(Death2010Onderzoeks)/nrow(onderzoeks)
verschilVerhoudingPCAdeathData <- nrow(PCAdeath2010Controle)/nrow(controle) nrow(PCAdeath2010Onderzoeks)/nrow(onderzoeks)
#### Analyse data dmv random groepen ####
toetsData <- toetsen(data2, 15871, 15867)
#### Plots van toets data ####
# Weergave van log-rank toets
plotToetsData(toetsData$logRank, LogRankData$chisq)
# Weergave van verschil percentage voor overleden mannen
plotToetsData(toetsData$verschilVerhoudingDeath, verschilVerhoudingDeathData)
# Weergave van verschil percentage voor overleden mannen aan prostaatkanker
plotToetsData(toetsData$verschilVerhoudingPCAdeath,
verschilVerhoudingPCAdeathData)
67
68
D.2
BIJLAGE D. ANALYSE MET VERHOUDINGEN
Gesimuleerde datasets
#### Analyse dataset 1 ####
source(file="functies/bep_functies.r")
load("Simulatie/dataSimulatieControle30000.Rdata")
#### Bepalen verhoudingen verschil ####
controleSimulatie <- subset(dataSimulatieControle, group == 0)
onderzoeksSimulatie <- subset(dataSimulatieControle, group == 1)
#Dataset met aantal personen dat overleden is per groep
Death2010ControleSimulatie <- subset(controleSimulatie, Death2010 == 1)
Death2010OnderzoeksSimulatie <- subset(onderzoeksSimulatie, Death2010 == 1)
# Bepaal verschil verhoudingen
verschilVerhoudingDeathSimulatie <- nrow(Death2010ControleSimulatie)/
nrow(controleSimulatie) nrow(Death2010OnderzoeksSimulatie)/
nrow(onderzoeksSimulatie)
#### Analyse dataset 1 dmv random groepen ####
toetsDataSimControle <- toetsen(dataSimulatieControle, 15000, 15000)
#### Plots van toets data van simulatie 1 met effect ####
# Weergave van log-rank toets
plotToetsData(toetsDataSimControle$logRank, LogRankSimulatie$chisq)
# Weergave van verschil percentage voor overleden personen
plotToetsData(toetsDataSimControle$verschilVerhoudingDeath,
verschilVerhoudingDeathSimulatie)
#### Analyse dataset 2 ####
load("Simulatie/dataSimulatie30000.Rdata")
#### Bepalen percentage verschil ####
controleSimulatie2 <- subset(dataSimulatie, group == 0)
onderzoeksSimulatie2 <- subset(dataSimulatie, group == 1)
#Dataset met aantal personen dat overleden is per groep
Death2010ControleSimulatie2 <- subset(controleSimulatie2, Death2010 == 1)
Death2010OnderzoeksSimulatie2 <- subset(onderzoeksSimulatie2, Death2010 == 1)
# Bepaal verschil percentages
verschilVerhoudingDeathSimulatie2 <- nrow(Death2010ControleSimulatie2)/
nrow(controleSimulatie2) nrow(Death2010OnderzoeksSimulatie2)/
nrow(onderzoeksSimulatie2)
D.2. GESIMULEERDE DATASETS
#### Analyse dataset 2 dmv random groepen ####
toetsDataSim <- toetsen(dataSimulatie, 15000, 15000)
#### Plots van toets data van simulatie 2 zonder effect ####
# Weergave van log-rank toets
plotToetsData(toetsDataSim$logRank, LogRankSimulatie2$chisq)
# Weergave van verschil percentage voor overleden mannen
plotToetsData(toetsDataSim$verschilVerhoudingDeath,
verschilVerhoudingDeathSimulatie2)
69
70
BIJLAGE D. ANALYSE MET VERHOUDINGEN
Bijlage E
Survivalanalyse
E.1
Prostaatkanker-dataset
load("Data/Dataset/data2.R")
#### Kaplan-Meier schatter
# Plotten van controle groep en onderzoeksgroep in n figuur (gehele groepen)
survivalTotaal <- survfit(Surv(ftime, Death2010)~group, data = data2,
conf.type="none")
plot(survivalTotaal, xlab="Tijd in jaren", ylab="Survival Kans", ymin = 0.5,
col=c("blue","red"))
legend(0.5,0.6, c("Controlegroep","Onderzoeksgroep"), lty=c(1,1),
lwd = c(2.5,2.5) ,col = c("blue","red"))
# Uitvoeren van de log-rank toets
LogRankData <- survdiff(Surv(ftime, Death2010)~group, data = data2, rho = 0)
#### Plots per leeftijd per groep met Kaplan-Meier schatter
#Plotten van controle groep en onderzoeksgroep in n figuur (gehele groepen)
kans = 0.90
tijdControles <- tijdenSurvivalKans(data2, 1, kans)
tijdInterventies <- tijdenSurvivalKans(data2, 2, kans)
# Plot maken van alle leeftijden tegen de tijdControles en tijdInterventies
plot(tijdControles, type = "l", xlim=c(54,75), ylim=c(0,15), col="blue")
par(new=T)
plot(tijdInterventies,type="l", xlim=c(54,75), ylim=c(0,15),col="red")
par(new=F)
legend(69,15, c("Controlegroep","Onderzoeksgroep"), lty=c(1,1),
lwd = c(2.5,2.5) ,col = c("blue","red"))
71
72
E.2
BIJLAGE E. SURVIVALANALYSE
Gesimuleerde datasets
source(file="functies/bep_functies.r")
#### Analyse dataset 1 ####
# Inlezen van de gesimuleerde data
load("Simulatie/dataSimulatieControle30000.Rdata")
#### Kaplan-Meier schatter
# Plotten van controle groep en onderzoeksgroep in n figuur (gehele groepen)
survivalTotaal <- survfit(Surv(ftime, Death2010)~group,
data = dataSimulatieControle,
conf.type="none")
plot(survivalTotaal, xlab="Tijd", ylab="Survival Kans", ymin = 0.5,
col=c("blue","red"))
legend(0.5,0.6, c("Controlegroep","Onderzoeksgroep"), lty=c(1,1),
lwd = c(2.5,2.5) ,col = c("blue","red"))
LogRankSimulatie <- survdiff(Surv(ftime, Death2010)~group,
data = dataSimulatieControle, rho = 0)
#### Plots per leeftijd per groep met Kaplan-Meier schatter
kans = 0.85
tijdControles <- tijdenSurvivalKans(dataSimulatieControle, 1, kans)
tijdInterventies <- tijdenSurvivalKans(dataSimulatieControle, 2, kans)
# Plot maken van alle leeftijden tegen de tijdControles en tijdInterventies
plot(tijdControles, type = "l", xlim=c(54,75), ylim=c(0,15), col="blue")
par(new=T)
plot(tijdInterventies,type="l", xlim=c(54,75), ylim=c(0,15),col="red")
par(new=F)
legend(69,15, c("Controlegroep","Onderzoeksgroep"), lty=c(1,1),
#### Analyse dataset 2 ####
# Inlezen van de gesimuleerde data
load("Simulatie/dataSimulatie30000.Rdata")
#### Kaplan-Meier schatter
# Plotten van controle groep en onderzoeksgroep in n figuur (gehele groepen)
survivalTotaal <- survfit(Surv(ftime, Death2010)~group, data = dataSimulatie,
conf.type="none")
plot(survivalTotaal, xlab="Tijd", ylab="Survival Kans", ymin = 0.5,
col=c("blue","red"))
legend(0.5,0.6, c("Controlegroep","Onderzoeksgroep"), lty=c(1,1),
lwd = c(2.5,2.5) ,col = c("blue","red"))
LogRankSimulatie2 <- survdiff(Surv(ftime, Death2010)~group,
data = dataSimulatie, rho = 0)
E.2. GESIMULEERDE DATASETS
#### Plots per leeftijd per groep met Kaplan-Meier schatter
kans = 0.85
tijdControles <- tijdenSurvivalKans(dataSimulatie, 1, kans)
tijdInterventies <- tijdenSurvivalKans(dataSimulatie, 2, kans)
# Plot maken van alle leeftijden tegen de tijdControles en tijdInterventies
plot(tijdControles, type = "l", xlim=c(54,75), ylim=c(0,15), col="blue")
par(new=T)
plot(tijdInterventies,type="l", xlim=c(54,75), ylim=c(0,15),col="red")
par(new=F)
legend(69,15, c("Controlegroep","Onderzoeksgroep"), lty=c(1,1),
lwd = c(2.5,2.5) ,col = c("blue","red"))
73
74
BIJLAGE E. SURVIVALANALYSE
Bijlage F
Competing risk analyse
F.1
Prostaatkanker-dataset
source(file="functies/bep_functies.r")
load("Data/Dataset/data2.R")
# Bepalen van de cumulatieve incidentiefuncties
cum <- cuminc(data2$ftime, data2$status, data2$group)
print(cum)
# Plot cumulatieve incidentiefuncties
plot(cum, ylim = c(0,0.5), xlab = "Tijd in jaren", ylab = "Hazard",
lty = c(1,1,2,2), col = c("blue","red","blue","red"))
legend(0,1,c("Controlegroep beoogde gebeurtenis",
"Onderzoeksgroep beoogde gebeurtenis",
"Controle groep competing risks",
"Onderzoeksgroep competing risks"),
lty=c(1,1,2,2) ,col = c("blue","red","blue","red"))
# Plot cumulatieve incidentiefuncties voor prostaatkanker
plot(cum[1], xlim = c(0,17), ylim = c(0,0.01), col = "blue")
par(new=T)
plot(cum[2], xlim = c(0,17), ylim = c(0,0.01), col = "red")
par(new=F)
legend(0,0.01, c("Controlegroep","Onderzoeksgroep"), lty=c(1,1),
lwd = c(2.5,2.5) ,col = c("blue","red"))
F.2
Gesimuleerde datasets
load("Simulatie/dataSimulatieControle30000.Rdata")
#### Dataset 1 ####
# Bepalen van de cumulatieve incidentiefuncties
simCum1 <- cuminc(dataSimulatieControle$ftime,
dataSimulatieControle$status,
dataSimulatieControle$group)
75
76
BIJLAGE F. COMPETING RISK ANALYSE
# Plot cumulatieve incidentiefuncties
plot(simCum1, ylim = c(0,0.5), xlab = "Tijd in jaren", ylab = "Hazard",
lty = c(1,1,2,2), col = c("blue","red","blue","red"))
legend(0,0.5,c("Controlegroep beoogde gebeurtenis",
"Onderzoeksgroep beoogde gebeurtenis",
"Controle groep competing risks",
"Onderzoeksgroep competing risks"),
lty=c(1,1,2,2) ,col = c("blue","red","blue","red"))
# Plot cumulatieve incidentiefuncties voor de beoogde gebeurtenis
plot(simCum1[1], xlim = c(0,17), ylim = c(0,0.1), col = "blue")
par(new=T)
plot(simCum1[2], xlim = c(0,17), ylim = c(0,0.1), col = "red")
par(new=F)
legend(0,0.1, c("Controlegroep","Onderzoeksgroep"), lty=c(1,1),
lwd = c(2.5,2.5) ,col = c("blue","red"))
#### Dataset 2 ####
load("Simulatie/dataSimulatie30000.Rdata")
# Bepalen van de cumulatieve incidentiefuncties
simCum2 <- cuminc(dataSimulatie$ftime, dataSimulatie$status,
dataSimulatie$group)
# Plot cumulatieve incidentiefuncties
plot(simCum2, ylim = c(0,0.5), xlab = "Tijd in jaren", ylab = "Hazard",
lty = c(1,1,2,2), col = c("blue","red","blue","red"))
legend(0,0.5,c("Controlegroep beoogde gebeurtenis",
"Onderzoeksgroep beoogde gebeurtenis",
"Controle groep competing risks",
"Onderzoeksgroep competing risks"),
lty=c(1,1,2,2) ,col = c("blue","red","blue","red"))
# Plot cumulatieve incidentiefuncties voor de beoogde gebeurtenis
plot(simCum2[1], xlim = c(0,17), ylim = c(0,0.1), col = "blue")
par(new=T)
plot(simCum2[2], xlim = c(0,17), ylim = c(0,0.1), col = "red")
par(new=F)
legend(0,0.1, c("Controlegroep","Onderzoeksgroep"), lty=c(1,1),
lwd = c(2.5,2.5) ,col = c("blue","red"))
Bibliografie
[1] M. T. Koller, H. Raatz, E. W. Steyerberg, and M. Wolbers, “Competing risks and the
clinical community: irrelevance or ignorance?,” Statistics in Medicine, vol. 31, no. 11-12,
pp. 1089–1097, 2012.
[2] H. Putter, M. Fiocco, and R. B. Geskus, “Tutorial in biostatistics: competing risks and
multi-state models,” Statistics in Medicine, vol. 26, no. 11, pp. 2389–2430, 2007.
[3] K. C. Cain, S. D. Harlow, R. J. Little, B. Nan, M. Yosef, J. R. Taffe, and M. R. Elliott,
“Bias due to left truncation and left censoring in longitudinal studies of developmental and
disease processes,” American Journal of Epidemiology, vol. 173, no. 9, pp. 1078–1084, 2011.
[4] G. Jongbloed, Bio- and Environmental Statistics. 2013.
[5] D. G. Kleinbaum and M. Klein, Survival analysis. Springer, 1996.
[6] H. Putter, Dynamic Prediction in Clinical Survival Analysis. CRC Press, 2011.
[7] J. P. Fine and R. J. Gray, “A proportional hazards model for the subdistribution of a
competing risk,” Journal of the American Statistical Association, vol. 94, no. 446, pp. 496–
509, 1999.
[8] B. Haller, G. Schmidt, and K. Ulm, “Applying competing risks regression models: an
overview,” Lifetime Data Analysis, vol. 19, no. 1, pp. 33–58, 2013.
[9] M. J. Roobol, R. Kranse, C. H. Bangma, A. G. van Leenders, B. G. Blijenberg, R. H.
van Schaik, W. J. Kirkels, S. J. Otto, T. H. van der Kwast, H. J. de Koning, and F. H.
Schr¨
oder, “Screening for prostate cancer: Results of the rotterdam section of the european
randomized study of screening for prostate cancer,” European Urology, vol. 64, no. 4, pp. 530
– 539, 2013.
[10] B. Haller and K. Ulm, “Flexible simulation of competing risks data following prespecified
subdistribution hazards,” Journal of Statistical Computation and Simulation, no. ahead-ofprint, pp. 1–20, 2013.
77

Adequaat voorschrijfgedrag bij ouderen. Implementatie van

Programma – Dag van de Vriendenkring van de Logistiek 2014

Korte Inhoud en CV

vragen - Hogeschool Rotterdam

Ervaring met de vacuümpomp

Voorbeeldtoets Wiskunde

wetenschap versus veldrealiteit Huisarts-in-opleiding:Els

schrijftherapie WriteJunior - Expertisecentrum Forensische Jeugdzorg

2 x 2 jaar - Universiteit Hasselt

PTS POBOS TRAUMA SCAN www.p o bos.be