VŠB - Technická univerzita Ostrava Fakulta elektrotechniky a informatiky Katedra aplikované matematiky
Diplomová práce
2014
Michal Běloch
VŠB - Technická univerzita Ostrava Fakulta elektrotechniky a informatiky Katedra aplikované matematiky
Vybrané statistické nástroje pro verifikaci logistických regresních modelů Advanced statistical tools for verification of the logistic regression models
2014
Michal Běloch
- 2 -
- 3 -
Zadání diplomové práce Bc. Michal Běloch
Student: Studijní program:
1103T031 Výpočetní matematika
Studijní obor: Téma:
N2647 Informační a komunikační technologie
Vybrané statistické nástroje pro verifikaci logistických regresních modelů. Advanced statistical tools for verification of the logistic regression models.
Zásady pro vypracování: Logistická regrese je zvláštní zejména tím, že vysvětlovaná proměnná je zde binární diskrétní proměnná. Nejčastěji je logistická regrese využívána pro účely zpracování lékařských dat. Cílem této práce je vytvoření matematického aparátu (zejména jeho PC implementace) pro posouzení vhodnosti buď předpokládaných nebo nově generovaných logistických regresních modelů. Postup práce: 1. Studium základů logistické regresní analýzy. 2. Předběžné posouzení logistického modelu pomocí lineární a exponenciální analýzy, algoritmizace. 3. Studium statistických testů pro posouzení vhodnosti modelů v logistické regresi. 4. Algoritmické zpracování a počítačová implementace vybraných testů této kategorie. 5. Provedení a vyhodnocení těchto testů na vzorku lékařských dat dle instrukcí vedoucího práce. Seznam doporučené odborné literatury: 1. Hosmer D.W., Lemeshow S., Applied Logistic Regression, Wiley 2000, ISBN 0-47135632-8. 2. Hebák P., Hustopecký J., Malá I., Vícerozměrné statistické metody (2), INFORATORIUM 2005, ISBN 80-7333-036-9.
- 4 -
Prohlašuji, že jsem tuto diplomovou práci vypracoval samostatně. Uvedl jsem všechny literární prameny a publikace, ze kterých jsem čerpal.
V Ostravě dne 07. května 2014
.............................................
- 5 -
Rád bych zde poděkoval svému vedoucímu prof. Ing. Radimu Brišovi, Csc. za poskytnutá data, literaturu a odbornou pomoc a rovněž své rodině a přátelům za snahu pochopit vše, co dělám.
- 6 -
Abstrakt První část této práce se zabývá tvorbou logistických regresních modelů za účelem zhodnocení morbidity pacientů s rakovinou kolorekta. Druhá část se zaměřuje na metody verifikace logistických regresních modelů a jejich algoritmizaci v programovacím jazyce R.
Klíčová slova exponenciální analýza, lineární analýza, logistická regrese, morbidita, R, verifikace modelu
Abstract First part of this work talks about logistic regression and its use in evaluation of morbidity of colorectal cancer patients. Second part focuses on verification methods for logistic regression and creation of algorithms in R programming language.
Keywords exponential analysis, linear analysis, logistic regression, morbidity, R, model verification
- 7 -
Seznam použitých symbolů a zkratek FNO MMV POSSUM
Fakultní nemocnice Ostrava – Poruba Metoda maximální věrohodnosti Physiological and Operative Severity Score for the enUmeration of Mortality and Morbidity
Seznam obrázků obr. 1: Výstup algoritmu pro exponenciální analýzu v R .................................................... 34 obr. 2: Koláčový graf pro rozložení pohlaví datového souboru .......................................... 36 obr. 3: Koláčový graf pro způsob operace ........................................................................... 37 obr. 4: Krabicový graf zachycující věkové rozložení pacientů ........................................... 37 obr. 5: Výskyt komplikací v závislosti na FS a OS u pacientů operovaných laparoskopicky ................................................................................................................................. 43 obr. 6: Výskyt komplikací v závislosti na FS a OS u pacientů operovaných otevřenou metodou ................................................................................................................... 43
- 9 -
Seznam tabulek tab. 1: Možnosti při testování hypotéz .............................................................................. 16 tab. 2: Rozhodování na základě p-value ........................................................................... 16 tab. 3: Tabulka lineární analýzy ........................................................................................ 24 tab. 4: Vzorek dat .............................................................................................................. 25 tab. 5: Příklad lineární analýzy pro vzorek dat ................................................................. 25 tab. 6: Příklad exponenciální analýzy pro vzorek dat ....................................................... 27 tab. 7: Příklad vylepšené exponenciální analýzy pro vzorek dat ...................................... 29 tab. 8: Rozložení veličin pro danou metodiku .................................................................. 38 tab. 9: Získání operačního skóre ....................................................................................... 38 tab. 10: Získání fyziologického skóre ............................................................................... 39 tab. 11: Lineární analýza pro pacienty operované laparoskopicky ................................... 40 tab. 12: Lineární analýza pro pacienty operované otevřenou metodou ............................ 40 tab. 13: Exponenciální analýza pro pacienty operované laparoskopicky ......................... 41 tab. 14: Exponenciální analýza pro pacienty operované otevřenou metodou ................... 41 tab. 15: P-value pro jednotlivé testy ................................................................................. 42 tab. 16: Lineární analýza pro pacienty operované laparoskopicky ................................... 45 tab. 17: Exponenciální analýza pro pacienty operované laparoskopicky ......................... 45 tab. 18: Vylepšená exponenciální analýza pro pacienty operované laparoskopicky ........ 46 tab. 19: P-value pro jednotlivé testy ................................................................................. 46 tab. 20: Lineární analýza pro pacienty operované otevřenou metodou ............................ 47 tab. 21: Exponenciální analýza pro pacienty operované otevřenou metodou ................... 47 tab. 22: Vylepšená exponenciální analýza pro pacienty operované otevřenou metodou .. 48 tab. 23: P-value pro jednotlivé testy ................................................................................. 48
- 10 -
Obsah 1.
Úvod ....................................................................................................................... 13
2.
Logistická regrese .................................................................................................. 14 2.1.
Testování hypotéz ........................................................................................... 14
2.1.1. Testové kritérium ....................................................................................... 15 2.1.2. Chyba I. a II. druhu .................................................................................... 15 2.1.3. Čistý test významnosti ............................................................................... 16 2.2.
Metoda maximální věrohodnosti .................................................................... 16
2.3.
Jednorozměrná logistická regrese ................................................................... 17
2.3.1. Testování významnosti koeficientů ............................................................ 19 2.4.
Vícerozměrná logistická regrese .................................................................... 21
2.4.1. Testování významnosti koeficientů ............................................................ 22 3.
4.
5.
Verifikace modelu .................................................................................................. 23 3.1.
Lineární analýza ............................................................................................. 23
3.2.
Exponenciální analýza .................................................................................... 25
3.3.
Vylepšená exponenciální analýza ................................................................... 28
3.4.
Pearsonova statistika a deviance..................................................................... 29
3.5.
Hosmer – Lemeshowy testy ........................................................................... 31
Popis algoritmů ...................................................................................................... 33 4.1.
R...................................................................................................................... 33
4.2.
Požadavky na vstup ........................................................................................ 33
4.3.
Lineární a exponenciální analýza ................................................................... 34
4.4.
Funkce vracející p-value ................................................................................. 35
Zpracování dodaných dat ....................................................................................... 36 5.1.
Datový soubor................................................................................................. 36
5.1.1. Operační a fyziologické skóre, model POSSUM ....................................... 38 5.2.
Aplikace známých modelů ............................................................................. 39
5.3.
Vytvoření nového modelu .............................................................................. 42
5.4.
Verifikace nového modelu.............................................................................. 44
5.4.1. Verifikace modelu pro laparoskopickou operaci ....................................... 45 5.4.2. Verifikace modelu pro otevřenou operaci .................................................. 47 - 11 -
6.
Závěr ...................................................................................................................... 49
7.
Literatura a reference ............................................................................................. 50
- 12 -
Úvod ___________________________________________________________________________
1. Úvod Rakovina kolon a rekta, obvykle souhrnně označována jako rakovina kolorekta, patří do pětice nejčastějších rakovinových onemocnění jak u nás, tak i ve světě. Z tohoto důvodu je tomuto typu rakoviny věnována vysoká pozornost většiny lékařských zařízení. Obrovský pokrok medicíny na přelomu tisíciletí vedl k vyvinutí léčebných postupů, které umožňují úspěšnou léčbu této zákeřné choroby. Jednou takovou metodou je laparoskopická operace. Tato moderní operační metoda si postupně získala značnou oblibu hlavně proto, že přináší nejmenší zátěž pro organismus pacienta. V případě operací nás v konečném důsledku zajímají dvě hodnoty: mortalita a morbidita. Z hlediska mortality je obvykle sledována časná mortalita, obvykle do třiceti dnů od provedení operace. Sledování morbidity je obvykle komplikovanější. Komplikací, které se mohou vyskytnout, je široké spektrum a ne vždy musí být přímo způsobeny chirurgickým zásahem. Zřejmě i z tohoto důvodu se procentuální výskyt komplikací u provedených operací vyskytuje v rozmezí od 4 % až po 26 %. Pro lékaře by byl velmi užitečný nástroj, který by byl schopen určit riziko výskytu komplikací. Mnoho lékařských pracovišť se touto problematikou zabývá a již došlo k vytvoření mnoha skórovacích systémů. Jedním z takovýchto systémů je POSSUM (Physiological and Operative Severity Score for the enUmeration of Mortality and Morbidity). Tento model byl vytvořen již v roce 1991 a je využíván nejen pro odhad rizik u operací kolorekta, ale i u jiných typů zákroků, a je hojně aplikován zejména ve Velké Británii.
- 13 -
Logistická regrese ___________________________________________________________________________
2. Logistická regrese Logistická regrese je dnes velice populární metoda. Přes skromné počátky v epidemiologickém výzkumu se postupně dostala do dalších odvětví a v současnosti se s ní můžeme setkat v biomedicíně, ekonomice, kriminalistice, ekologii a sociologii. V této kapitole nadefinujeme jednorozměrnou logistickou regresi, kterou pak dále rozšíříme na vícerozměrnou verzi. Rovněž se zde zaměříme na testování významnosti jednotlivých koeficientů. Protože budeme dále používat testování hypotéz a metodu maximální věrohodnosti, nadefinujeme nejdříve tyto pojmy.
2.1. Testování hypotéz Statistickou hypotézou rozumíme výrok o rozdělení pozorované náhodné veličiny. Tato hypotéza může pojednávat o parametrech rozdělení náhodné veličiny nebo o vlastnostech náhodné veličiny. Test statistické hypotézy implikuje rozhodovací proces, kdy na základě výběrového souboru rozhodujeme o nezamítnutí nulové hypotézy, nebo o zamítnutí nulové hypotézy ve prospěch alternativní hypotézy. Nulová hypotéza H0 vyjadřuje nulovost sledovaného efektu. Obvykle je vyjádřena rovností mezi testovaným parametrem θ a jeho očekávanou hodnotou θ0 . H 0 : θ = θ0 Alternativní hypotéza HA je námi vybraná tak, aby popírala tvrzení dané nulové hypotézy. To nám tedy dává celkem čtyři možnosti jak zformulovat alternativní hypotézu.
1) H A : θ = θ1 2 ) H A : θ ≠ θ0 3) H A : θ < θ 0 4 ) H A : θ > θ0
- 14 -
Logistická regrese ___________________________________________________________________________ 2.1.1. Testové kritérium Obor hodnot testovaného parametru θ se dělí na dvě disjunktní množiny. První je obor přijetí V pro testovanou hypotézu H0 , druhý se nazývá kritický obor W pro zamítnutí hypotézy H0 . Tento kritický obor stanovujeme tak, aby pravděpodobnost výskytu pozorované hodnoty testovaného parametru θ v něm byla velmi malá. Hranici mezi kritickým oborem a oborem přijetí nazýváme kritická hodnota testu tkrit . Padne-li pozorovaná hodnota testovaného parametru θ do oboru přijetí V, nulovou hypotézu H0 nezamítáme. Padne-li do kritického oboru W, pak nulovou hypotézu H0 zamítáme ve prospěch alternativy HA . Kritický obor W lze popsat pomocí kritického oboru W* testového kritéria T ( X ) . Testové kritérium T ( X ) je výběrová charakteristika, která má vztah k nulové hypotéze H0. Za předpokladu platnosti nulové hypotézy H0 známe rozdělení tohoto kritéria. Následně padne-li pozorovaná hodnota testového kritéria T ( X ) do kritického oboru W*, zamítáme nulovou hypotézu H0. V opačném případě H0 nezamítáme. 2.1.2. Chyba I. a II. druhu Při testování hypotéz mohou nastat čtyři různé případy. Ke správnému rozhodnutí dospějeme tehdy, platí-li nulová hypotéza a my ji nezamítáme, nebo platí alternativní hypotéza a my nulovou hypotézu zamítneme. V prvním případě hovoříme o spolehlivosti testu a značíme ji 1 − α . Výhodou je, že si tuto hladinu určujeme sami, a je nám tedy předem známa. Hodnotu α obvykle volíme 0,05, což nám dává spolehlivost 95%. Ve druhém případě hovoříme o síle testu, značíme ji 1 − β . Chyb se dopouštíme v případě zamítnutí nulové hypotézy i tehdy, je-li ve skutečnosti platná. Tuto chybu nazýváme chyba I. druhu, dopouštíme se jí s pravděpodobností α a nazýváme ji hladina významnosti. Chybou II. druhu značíme rozhodnutí, kdy nulovou hypotézu nezamítáme i přesto, že ve skutečnosti platí alternativa. Této chyby se dopouštíme s pravděpodobností β . Žádná z těchto chyb nelze zcela eliminovat, proto se snažíme postupovat tak, abychom se dopouštěli co nejmenší chyby. Bohužel pokud snížíme β, zvyšujeme zároveň hladinu významnosti a naopak. Musíme tedy najít ideální poměr těchto chyb. Spolehlivost testu si volíme sami. Sílu testu můžeme zvýšit volbou vhodné testovací metody. Nejlepším způsobem, jak snížit β je ale zvýšení rozsahu našeho výběrového souboru. V tomto případě snižujeme chybu II. druhu, aniž bychom zvyšovali chybu I. druhu.
- 15 -
Logistická regrese ___________________________________________________________________________
Skutečnost
Naše zjištění Nezamítáme H0 Zamítáme H0 Platí H0
Neplatí H0
Správně 1-α: spolehlivost testu Chyba II. druhu β
Chyba I. druhu α: hladina významnosti Správně 1-β: síla testu
tab. 1: Možnosti při testování hypotéz
2.1.3. Čistý test významnosti Při testování hypotéz můžeme využít dvou přístupů. Prvním je klasický test, který je však v moderních aplikacích nahrazen druhým přístupem, čistým testem významnosti. U klasického testu je součástí vstupu i hladina významnosti α , zatímco u čistého testu tuto hodnotu a priori znát nemusíme. Postup při čistém testu významnosti je: 1) Formulace H0 a HA 2) Volba testového kritéria T ( X ) 3) Výpočet pozorované hodnoty xobs testové statistiky T ( X ) 4) Výpočet p-value 5) Rozhodnutí na základě p-value P-value je číselná hodnota, která nám říká, jaká je nejnižší hladina významnosti, na které můžeme nulovou hypotézu H0 zamítnout, a zároveň nejvyšší hladina významnosti, na níž se už nulová hypotéza H0 nezamítá. S klesající hodnotou p-value roste výpověď náhodného výběru proti nulové hypotéze H0. K rozhodování tedy využíváme následující schéma: P-value <α >α
Rozhodnutí Zamítnutí H0 ve prospěch HA Nezamítnutí H0
tab. 2: Rozhodování na zá kladě p-value
2.2. Metoda maximální věrohodnosti MMV je jednoduchá metoda používaná pro konstrukci odhadů neznámých parametrů známých rozdělení pravděpodobnosti. Při MMV maximalizujeme věrohodnostní funkci, která představuje sdruženou hustotu pravděpodobnosti daného náhodného výběru, chápanou jako funkci neznámých parametrů. Odhady získané touto metodou se obvykle vyznačují dobrými statistickými vlastnostmi. - 16 -
Logistická regrese ___________________________________________________________________________
T
Nechť ( t1 ,..., tn ) je náhodný výběr z rozdělení s hustotou f ( t; Θ ) , kde Θ je neznámý parametr. Snažíme se najít tzv. věrohodnostní funkci danou n
L ( t1 ,..., tn ; Θ ) = f ( t1 ; Θ ) ⋅ f ( t2 ; Θ ) ⋅ ... ⋅ f ( tn ; Θ ) = ∏ f ( ti ; Θ ) i =1
ˆ tak, aby Θ ˆ =Θ ˆ ( t ,..., t ) bylo co nejlepším odhadem pro Θ . Pravá strana a z ní získat Θ 1 n rovnice je sdružená hustota pravděpodobnosti n-nezávislých proměnných ( t1 ,..., tn ) se stejným rozdělením. Věrohodnostní funkce L je funkcí neznámého parametru Θ , který je odhadován. Metoda maximální věrohodnosti je založena na nalezení takové hodnoty Θ , aby hodnota věrohodnostní funkce L byla maximální. Praktické aplikace ukázaly, že je výhodnější maximalizovat funkci ln L . Maximálně věrohodným odhadem parametru Θ nazveme hodnotu parametru získanou ∂ ln L ( t1 ,..., tn ; Θ ) z této rovnice: = 0. ∂Θ Pokud je hledaných parametrů více, přejde výše uvedená rovnice na soustavu rovnic ∂ ln L ( t1 ,..., tn ; Θ1 ,..., Θ p ) ∂Θi
= 0 pro i = 1, 2,..., p a jejím vyřešením získáme hledané parametry
Θ1 , Θ 2 ,..., Θ p .
2.3. Jednorozměrná logistická regrese Cílem modelování je vždy nalézt model, který nejlépe popisuje zdrojová data a jejich vztah ke konkrétní vysvětlované veličině. V logistické regresi je vysvětlovaná (závislá, outcome) veličina binární, nejčastěji nabývající hodnot 0 a 1, jev nenastal, nebo nastal. Mějme vzorek n nezávislých dvojic ( xi , yi ) , i = 1, 2,..., n , kde yi označuje binární závislou proměnnou a xi je hodnota nezávislé proměnné pro i-tý objekt. Model pro logistickou regresi vypadá následovně:
e β0 + β1x π ( x) = . 1 + e β0 + β1x V logistické regresi navíc využíváme tzv. logitovou transformaci, která je definována:
- 17 -
Logistická regrese ___________________________________________________________________________
π ( x) g ( x ) = ln = β 0 + β1 x . 1 − π x ( ) Tato transformace má řadu příjemných vlastností, které sdílí s lineární regresí. Logit
g ( x ) je lineární ve svých parametrech, může být spojitý a v závislosti na rozsahu x může nabývat hodnot od −∞ až po +∞ . Pro naplnění logistického modelu potřebujeme znát hodnoty parametrů β 0 , β1 . Pro jejich nalezení se v případě logistické regrese používá metoda maximální věrohodnosti. Pokud Y nabývá hodnot 0 a 1, dává výraz pro π ( x ) podmíněnou pravděpodobnost
P (Y = 1| x ) . Odtud snadno plyne, že výraz 1 − π ( x ) dává pravděpodobnost P (Y = 0 | x ) . Tedy pro dvojice ( xi , yi ) , kde yi = 1 , je přínos k věrohodnostní funkci π ( xi ) a pro dvojice, kde yi = 0 , je přínos 1 − π ( xi ) . Tento fakt můžeme jednoduše zapsat do vztahu 1− yi
yi
π ( xi ) 1 − π ( xi )
.
Protože předpokládáme nezávislost jednotlivých dvojic, dostaneme finální věrohodnostní funkci jako součin výše uvedeného výrazu: n
1− yi
yi
l ( β 0 , β1 ) = ∏ π ( xi ) 1 − π ( xi )
.
i =1
Pro maximalizaci použijeme funkci: n
yi
L ( β 0 , β1 ) = ln l ( β 0 , β1 ) = ln ∏ π ( xi ) 1 − π ( xi )
1− yi
=
i =1
n n e β0 + β1xi = ∑ yi ln π ( xi ) + (1 − yi ) ln 1 − π ( xi ) = ∑ yi ln β0 + β1 xi i =1 i =1 1+ e
n
e β0 + β1 xi y + 1 − ln 1 − ( ) i β0 + β1 xi 1+ e
= ∑ yi β 0 + β1 xi − ln 1 + e β0 + β1xi − (1 − yi ) ln 1 + e β0 + β1xi = y =1
(
)
(
n
)
(
= ∑ yi ( β 0 + β1 xi ) − ln 1 + e β0 + β1xi i =1
- 18 -
).
=
Logistická regrese ___________________________________________________________________________ Posledním krokem pro získání maximálních věrohodných odhadů pomocí MMV je sestavení věrohodnostních rovnic a jejich vyřešení. K tomu potřebujeme příslušné parciální derivace podle hledaných parametrů β 0 , β1 :
∂L ( β 0 , β1 ) ∂β 0 ∂L ( β0 , β1 ) ∂β1
e β0 + β1xi = ∑ yi − 1 + e β0 + β1xi i =1 n
n = ∑ yi − π ( xi ) i =1
xi eβ0 + β1xi n = ∑ xi yi − = ∑ xi yi − π ( xi ) . 1 + e β0 + β1xi i =1 i =1 n
Příslušná soustava věrohodnostních rovnic je: n
∑ y − π ( x ) = 0 i
i
i =1 n
∑ x y − π ( x ) = 0 i
i
.
i
i =1
Obdrželi jsme nelineární soustavu dvou rovnic o dvou neznámých. K vyřešení této soustavy se obvykle používají iterační metody, které jsou zpravidla již implementovány v příslušných statistických balíčcích. My se těmito metodami zabývat nebudeme. Po vyřešení této soustavy dostaneme maximálně věrohodný odhad β = ( β 0 , β1 ) , který budeme značit βˆ . Obdobně maximálně věrohodný odhad π ( xi ) bude značen πˆ ( xi ) .
2.3.1. Testování významnosti koeficientů Jakmile obdržíme odhad koeficientů, zaměří se naše pozornost na jednotlivé proměnné v modelu. To obvykle zahrnuje formulaci a testování statistické hypotézy za účelem zjištění, zda-li nezávislé proměnné v modelu mají významný vliv na vysvětlovanou proměnnou. Jeden z přístupů k testovaní významnosti koeficientů klade otázku, řekne-li nám model, který obsahuje danou proměnnou, o vysvětlované proměnné více než model, který tuto proměnnou neobsahuje. Tuto otázku zodpovíme porovnáním pozorovaných hodnot závislé proměnné s dvěma predikovanými hodnotami z modelu s a bez proměnné, která nás zajímá. Pokud je předpovídaná proměnná v jistém smyslu lepší nebo přesnější v modelu obsahujícím zkoumanou proměnnou než v modelu, který tuto proměnnou neobsahuje, považujeme tuto proměnnou za významnou. - 19 -
Logistická regrese ___________________________________________________________________________
Metodikou tedy bude porovnání pozorovaných hodnot závislé proměnné s predikovanými hodnotami, které získáme z modelu obsahující zkoumanou proměnnou, s modelem, který tuto proměnnou neobsahuje. V logistické regresi je toto porovnání založeno na logaritmické věrohodnostní funkci L ( β0 , β1 ) . Konceptuálně uvažujeme o pozorovaných hodnotách závislé proměnné jako o predikovaných hodnotách saturovaného modelu. Saturovaným modelem pak rozumíme model, který obsahuje tolik parametrů β , kolik je v souboru pozorování. Toto porovnání pak zachycuje tento vztah: věrohodnost naplňovaného modelu D = −2 ln věrohodnost saturovaného modelu
Hodnota uvnitř závorek se nazývá věrohodnostní poměr. Použitím záporného dvojnásobku logaritmu obdržíme hodnotu se známým rozdělením vhodným pro testování hypotéz. Dosazením obdržíme:
πˆ ( xi ) 1 − πˆ ( xi ) D = −2∑ yi ln + (1 − yi ) ln 1 − yi yi i =1 n
Statistice D se říká deviance. Navíc věrohodnost saturovaného modelu je rovna 1, což přímo plyne z definice saturovaného modelu, kde πˆ ( xi ) = yi a pro věrohodnost platí n
1− yi
∏ y × (1 − y ) yi i
i
= 1 . Vztah pro devianci D tedy můžeme zjednodušit:
i =1
D = −2 ln [ věrohodnost naplňovaného modelu ] . Pro zhodnocení významnosti nezávislé proměnné porovnáme hodnotu deviance D s a bez nezávislé proměnné v rovnici. Změnu v D pak vyjádříme:
G = D [ model bez proměnné] − D [ model s proměnnou ] . Věrohodnost saturovaného modelu je obsažena v obou hodnotách D ve zmíněném rozdílu, proto můžeme vztah pro G zavést následovně:
věrohodnost bez proměnné G = −2 ln . věrohodnost s proměnnou
- 20 -
Logistická regrese ___________________________________________________________________________ Pro speciální případ jedné nezávislé proměnné se dá ukázat, že v případě absence n proměnné v modelu je maximálně věrohodný odhad β 0 roven ln 1 , kde n1 = ∑ yi , n0
n0 = ∑ (1 − yi ) a predikovaná hodnota je konstantní
n1 . V tomto případě pro G platí: n
n0 n1 n1 n0 n n G = −2 ln n yi 1− yi ∏ πˆ ( xi ) 1 − πˆ ( xi ) i =1
=
n
= 2∑ yi ln πˆ ( xi ) + (1 − yi ) ln 1 − πˆ ( xi ) − 2 [ n1 ln n1 + n0 ln n0 − n ln n ] i =1
Pokud testujeme hypotézu β1 = 0 , pak G ∼ χ 2 (1) ( G má chi-kvadrát rozdělení s jedním stupněm volnosti). Pokud používáme nějaký statistický software, je v něm obvykle jistá forma testu významnosti koeficientů implementována a jako výstup obdržíme p-value. V případě programu Statgraphics testujeme hypotézu:
H0: Výpověď testovaného parametru není dostatečně významná a daný parametr můžeme z modelu vynechat HA: Výpověď testovaného parametru je statisticky významná a měl by být ponechán Rozhodnutí pak učiníme na základě obdržené p-value v souladu s postupem při čistém testu významnosti.
2.4. Vícerozměrná logistická regrese Mějme p nezávislých proměnných vektorově zapsaných x = ( x1 , x2 ,..., x p ) . Označme podmíněnou pravděpodobnost, že výstup nastal P (Y = 1| x ) = π ( x ) . Logit vícerozměrného logistického regresního modelu má tvar:
g ( x ) = β 0 + β1 x1 + β 2 x2 + ... + β p x p . Příslušný model logistické regrese je pak ve tvaru: - 21 -
Logistická regrese ___________________________________________________________________________
e() π ( x) = g x . 1+ e ( ) g x
Předpokládejme, že máme vzorek n nezávislých pozorování ( xi , yi ) , i = 1, 2,..., n . Stejně jako v jednorozměrném případě potřebujeme pro naplnění regresního modelu získat odhady vektoru β = ( β 0 , β1 , β 2 ,..., β p ) . I ve vícerozměrném případě k jejich získání využijeme metodu maximální věrohodnosti. Věrohodnostní funkce bude v podstatě identická jako v jednorozměrném případě, jediným rozdílem je nepatrně rozdílná definice π ( x ) . Po derivování log věrohodnostní funkce vzhledem k p + 1 parametrům obdržíme soustavu p + 1 věrohodnostních rovnic, která vypadá následovně: n
∑ y − π ( x ) = 0 i
i
i =1 n
∑x
ij
yi − π ( xi ) = 0, j = 1, 2,..., p
i =1
Řešení této soustavy je opět numerické povahy a příslušná metodologie je již zabudována ve statistických balíčcích. Označme řešení této soustavy βˆ , získané hodnoty z modelu pro logistickou regresi πˆ ( xi ) .
2.4.1. Testování významnosti koeficientů Stejně jako v jednorozměrné logistické regresi je prvním krokem po naplnění modelu zvážení významnosti jednotlivých proměnných. Test věrohodnostního poměru pro celkovou významnost p koeficientů pro nezávislé proměnné modelu je prováděn naprosto stejně a je založen na hodnotě G jako v jednorozměrném případě. V tomto případě hledané hodnoty πˆ obsahují vektor s p + 1 parametry βˆ . Stanovíme-li nulovou hypotézu tak, že p hodnot koeficientů „zakřivení“ v modelu jsou nulové, má G chi-kvadrát rozdělení s p stupni volnosti: G ∼ χ 2 ( p ) .
- 22 -
Verifikace modelu ___________________________________________________________________________
3. Verifikace modelu V předchozí kapitole jsme vytvořili regresní model a rozhodli jsme o tom, které proměnné jsou statisticky významné. V této části se zaměříme, jak dobře model popisuje vysvětlovanou proměnnou. Nejdříve se podíváme na lineární a exponenciální analýzu. Tyto metody nám dávají zběžný pohled na to, jak dobře model predikuje, aniž bychom se museli zabývat testováním hypotéz. Použití těchto metod bylo navrženo v [9]. Druhou skupinou budou verifikační metody, kde již přikročíme k testování hypotéz, a tedy výpověď těchto metod bude silnější.
eg( x) máme V této kapitole předpokládáme, že logistický regresní model π ( x ) = 1 + e g ( x) k dispozici, tedy pro logitovou funkci g ( x ) = β 0 + β1 x1 + ... + β p x p známe hodnotu všech parametrů β 0 , β1 ,..., β p .
3.1. Lineární analýza Lineární analýza je jednoduchá metoda, jejíž výstup se dá nejjednodušeji charakterizovat následující tabulkou: Procentuální skupina [%]
Počet Počet výskytů pozorování sledovaného ve skupině znaku ve skupině
Předpovídaný počet výskytů sledovaného znaku ve skupině
<10
n1
v1
p1 = 0, 045* n1
10 – 19
n2
v2
p2 = 0,145* n2
20 – 29
n3
v3
p3 = 0, 245* n3
30 – 39
n4
v4
p4 = 0,345* n4
40 – 49
n5
v5
p5 = 0, 445* n5
- 23 -
Poměr sledovaných znaků ku předpovídaným znakům v1 p1 v2 p2 v3 p3 v4 p4 v5 p5
Verifikace modelu ___________________________________________________________________________ 50 – 59
n6
v6
p6 = 0,545* n6
60 – 69
n7
v7
p7 = 0, 645* n7
70 – 79
n8
v8
p8 = 0, 745* n8
80 – 89
n9
v9
p9 = 0,845* n9
>89
n10
v10
p10 = 0,945* n10
Celkem
n = ∑ ni
10
i =1
10
v6 p6 v7 p7 v8 p8 v9 p9 v10 p10
10
v = ∑ vi
v p
p = ∑ pi
i =1
i =1
tab. 3: Tabul ka lineární analýzy
První sloupec je pevně dán. Hodnoty ve druhém sloupci obdržíme tak, že dosadíme do
eg( x) modelu π ( x ) = hodnoty vysvětlujících proměnných x = ( x1 , x2 ,..., x p ) a podle 1 + e g ( x) získaného výsledku zvedneme počítadlo n v příslušném řádku. Pokud navíc tato konkrétní realizace náhodného vektoru měla přítomný znak vysvětlované proměnné, zvedneme počítadlo i ve třetím sloupci. V pátém sloupci nás zajímá poměr počtu skutečných výskytů a počtu predikovaného výskytu. Ideálně bychom si přáli, aby na každém řádku bylo číslo blízké jedničce, tedy skutečnost se shoduje s předpovědí. Číslo větší než jedna značí optimistický model, kdy výskyt vysvětlované proměnné je ve skutečnosti vyšší, než náš model předpovídá. Hodnoty menší než jedna pak značí pesimistický model, kdy model předpovídá větší výskyt vysvětlované proměnné, než je skutečnost (předpokládáme, že přítomnost znaku u vysvětlované proměnné je negativum, např. morbidita, mortalita atd.). Příklad: Mějme danou logit funkci g ( x ) = −1, 434 + 0, 025* FS + 0, 054* OS a vzorek dat s vysvětlovanou proměnnou Pooperační komplikace: Fyziologické Operační skore FS skore OS 12 0 15 11 15 20 18 10 20 10 11 17 23 8 15 13
Pooperační komplikace 1 0 1 0 0 0 1 0 - 24 -
π ( x ) *100 24 39 50 39 40 44 39 41
Verifikace modelu ___________________________________________________________________________ 19 14
17 9
0 0
49 35
tab. 4: Vzore k dat
Tabulka pro lineární analýzu bude vypadat následovně:
Skupina [%] <10 10 – 19 20 – 29 30 – 39 40 – 49 50 – 59 60 – 69 70 – 79 80 – 89 >89 Celkem
Predikovaný Počet Počet počet pozorování komplikací komplikací 0 0 0 0 0 0 1 1 0 4 1 1 4 0 2 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 10 3 4
Poměr komplikací ku předpovídaný 1 1 N 1 0 1 1 1 1 1 0,75
tab. 5: Pří klad lineární analýzy pro vzorek dat
Z příkladu vidíme, že v lineární analýze může nastat problém pro dělení nulou. V tomto 0 číslo ≠ 0 případě budeme brát výraz = 1 a =N. 0 0
3.2. Exponenciální analýza Podobně jako u lineární analýzy se u exponenciální analýzy snažíme vytvořit procentuální skupiny a porovnat počty výskytů vysvětlované proměnné s predikovanými hodnotami. Jednotlivá pozorování rozdělíme do skupin podle vypočtené hodnoty π ( x ) na základě jednotlivých realizací vysvětlujících proměnných, vyjádřených v procentech. Analýza se provádí od skupiny 90 – 100 a postupně zvětšujeme interval po 10 procentech (80 – 100, 70 – 100 atd.) až do doby, kdy bude vyhovovat podmínce, aby předpovídaný počet pacientů měl vzrůstající tendenci (byl stejný) vzhledem k předchozí skupině. Analýza je prováděna zdola nahoru. Jestliže dojde k porušení podmínky, analýza zdola nahoru se zastaví a skupina, u které došlo k porušení podmínky, se do analýzy nepočítá. Pokračuje se dále analýzou shora dolů. Pro tuto analýzu je důležitá poslední skupina analýzy zdola dolů, u které nedošlo
- 25 -
Verifikace modelu ___________________________________________________________________________ k porušení podmínky. Její dolní hranice bude určovat horní hranici skupin u analýzy shora dolů. Předpovídaný počet dostaneme následujícím způsobem: předpověď = počet ve skupině*dolní hranice skupiny/100
Příklad: Použijeme data z předchozího příkladu. Začínáme shora:
Skupina [%] 90 – 100
Počet Počet pozorování komplikací 0
0
Predikovaný počet komplikací 0
Rozšíříme interval:
Skupina [%] 80 – 100 90 – 100
Počet Počet pozorování komplikací 0 0
0 0
Predikovaný počet komplikací 0 0
Takto pokračujeme dále až obdržíme:
Skupina [%] 20 – 100 30 – 100 40 – 100 50 – 100 60 – 100 70 – 100 80 – 100 90 – 100
Počet Počet pozorování komplikací 10 9 5 1 0 0 0 0
3 2 1 1 0 0 0 0
Predikovaný počet komplikací 2 3 2 1 0 0 0 0
Vidíme, že u intervalu 20 – 100 došlo k porušení podmínky. Predikovaný počet komplikací je 2, avšak v předchozí skupině 30 – 100 je predikovaný počet komplikací 3. Interval 20 – 100 tedy do analýzy nebudeme počítat a začneme s novou horní mezí intervalu 30:
Skupina [%] 20 – 30 30 – 100
Počet Počet pozorování komplikací 1 9
1 2 - 26 -
Predikovaný počet komplikací 0 3
Verifikace modelu ___________________________________________________________________________ 40 – 100 50 – 100 60 – 100 70 – 100 80 – 100 90 – 100
5 1 0 0 0 0
1 1 0 0 0 0
2 1 0 0 0 0
Opět pokračujeme stejným způsobem, až obdržíme:
Skupina [%] 0 – 30 10 – 30 20 – 30 30 – 100 40 – 100 50 – 100 60 – 100 70 – 100 80 – 100 90 – 100
Počet Počet pozorování komplikací 1 1 1 9 5 1 0 0 0 0
1 1 1 2 1 1 0 0 0 0
Predikovaný počet komplikací 0 0 0 3 2 1 0 0 0 0
Finální tabulka exponenciální analýzy vypadá následovně:
Skupina [%] 0 – 30 10 – 30 20 – 30 30 – 100 40 – 100 50 – 100 60 – 100 70 – 100 80 – 100 90 – 100 Celkem
Predikovaný Počet Počet počet pozorování komplikací komplikací 1 1 0 1 1 0 1 1 0 9 2 3 5 1 2 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 10 3 3
Poměr komplikací ku předpovídaný N N N 0,66 0,5 1 1 1 1 1 1
tab. 6: Pří klad expone nciální analýzy pro vzorek dat
Nepraktickou vlastností exponenciální analýzy je nutnost několikerého přepočtu v případě porušení podmínky. Navíc se tento bod přepočtu liší soubor od souboru. Navíc je nejasné, jak přistupovat k analýze jednotlivých podskupin. Pokud bychom například chtěli zjistit poměr ve skupině 20 – 40, nemůžeme jednoduše sečíst hodnoty ve skupinách 20 – 30 a 30 – 40, protože zde došlo k přepočtu a pásmo 30 – 40 rovněž nemůžeme nezávisle - 27 -
Verifikace modelu ___________________________________________________________________________ analyzovat vůči pásmu 30 – 100. Druhým problémem je nutnost počítání a přepočítávání jednotek v daném pásmu. Jednotka, která je v oblasti 90 – 100 je započítána i pro 80 – 100, 70 – 100 až po 30 – 100, zatímco jiná jednotka s pravděpodobností 35% je započtena jednou jen v pásmu 30 – 100. I přes tyto nedostatky a přes velmi omezený datový soubor vidíme, že v celkovém součtu exponenciální analýza predikovala trošku lépe, než lineární.
3.3. Vylepšená exponenciální analýza U exponenciální analýzy jsme predikované hodnoty dostali takto: předpověď = počet ve skupině*dolní hranice skupiny/100. Jedná se však o odhad, který nemusí dostatečně popisovat chování modelu. Alternativou je použít větu o úplné pravděpodobnosti pro výpočet predikované hodnoty. Označme x modelem predikovaný počet komplikací ve skupině ( a, b ) , která obsahuje n pacientů a P ( A) pravděpodobnost, že náhodně vybraný pacient ze skupiny ( a, b ) bude mít komplikace. Potom x = nP ( A) . Jestliže se ve skupině ( a, b ) vyskytlo r různých pravděpodobností komplikací p1 , p2 ,..., pr a ni je počet pacientů s pravděpodobností n r r komplikace pi , pak můžeme psát x = n ∑ pi i = ∑ pi ni . Pokud provedeme přeznačení a i =1 n i =1 každému i-tému pacientovi ze skupiny ( a, b ) přiřadíme pravděpodobnost komplikace pi , pak dostaneme: n
x = ∑ pi i =1
Počet predikovaných komplikací v dané skupině pacientů tedy určíme jako součet pravděpodobností komplikací jednotlivých pacientů této skupiny.
Příklad: Na předchozí příklad exponenciální analýzy aplikujeme vylepšenou exponenciální analýzu:
Skupina [%] 0 – 30 10 – 30 20 – 30 30 – 100 40 – 100
Predikovaný Počet Počet počet pozorování komplikací komplikací 1 1 0 1 1 0 1 1 0,24 (0,2) 9 2 3,76 (2,7) 5 1 2,24 (2) - 28 -
Poměr komplikací ku předpovídaný N N N 0,53 0,44
Verifikace modelu ___________________________________________________________________________ 50 – 100 60 – 100 70 – 100 80 – 100 90 – 100 Celkem
1 0 0 0 0 10
1 0 0 0 0 3
1 (1) 0 0 0 0 3
1 1 1 1 1 1
tab. 7: Pří klad vylepšené exp onenciální analýzy pro vzorek dat
V tabulce je uveden předpovídaný počet komplikací pomocí věty o úplné pravděpodobnosti, v závorce jsou hodnoty z původní analýzy. I když na tomto příkladu vypadá vylepšená exponenciální analýza jako horší, musíme brát na vědomí extrémně nízký počet použitých dat. Data nebyla vybrána pro rozhodnutí nejvhodnějšího postupu, ale pro ilustraci vytváření jednotlivých analýz.
3.4. Pearsonova statistika a deviance V textu se dále objeví pojem kovariační schéma. Tímto rozumíme jednoznačnou kombinaci hodnot všech vysvětlujích proměnných. Máme-li dvě vysvětlující proměnné
{( 0, 0) , ( 0,1) , (1, 0) , (1,1)} . x = ( x , x ,..., x ) . Pokud
kódované {0,1} a {0,1} ,bude kovariační schéma vypadat takto: Dále označíme J jako počet unikátních pozorovaných hodnot
1
2
p
některé subjekty mají stejné hodnoty x , platí J < n . Počet subjektů s x = x j označíme jako J
m j pro j = 1, 2,..., J . Je zřejmé, že
∑m
j
= n . Obvykle předpokládáme, že J ≈ n , jelikož
j =1
očekáváme přítomnost alespoň jedné spojité vysvětlující proměnné v modelu. V logistické regresi je více způsobů, jak měřit rozdíl mezi pozorovanými a vypočtenými hodnotami. Pro zdůraznění, že vypočtené hodnoty jsou v logistické regresi počítány pro každé kovariační schéma, označíme j-té kovariační schéma jako yˆ j , kde
yˆ j = m jπˆ j = m j
e
( )
gˆ x j
1+ e
( )
gˆ x j
a gˆ ( x j ) rozumíme odhadnutý logit. Pro konkrétní kovariační schéma je Pearsonovo residuum definováno jako:
- 29 -
Verifikace modelu ___________________________________________________________________________
r ( y j , πˆ j ) =
(y
j
− m jπˆ j )
m jπˆ j (1 − πˆ j )
.
Statistika založená na těchto residuích se nazývá Pearsonova chí-kvadrát statistika: J
2
Χ 2 = ∑ r ( y j , πˆ j ) . j =1
Provádíme tedy součet residuí přes všechna kovariační schémata, která se nám vyskytují v datech. Residua deviance definujeme:
y mj − yj d ( y j , πˆ j ) = ± 2 y j ln j + ( m j − y j ) ln , m jπˆ j m j (1 − πˆ j ) znaménko před odmocninou je shodné se znaménkem výrazu ( y j − m jπˆ j ) . Pro kovariační schéma, kde by nastal případ, že y j = 0 :
d ( y j , πˆ j ) = − 2m j ln (1 − πˆ j ) a pro případ, že y j = m j :
d ( y j , πˆ j ) = 2m j ln (πˆ j ) . Součtová statistika založena na residuích deviancí je deviance J
D = ∑ d ( y j , πˆ j )
2
.
j =1
V případě, že J = n , se jedná o stejnou hodnotu, jaká byla vyjádřena pro devianci v kapitole 3.3.1.
- 30 -
Verifikace modelu ___________________________________________________________________________ Rozdělení statistiky Χ 2 a D za předpokladu, že získaný model je korektní ve všech předpokladech, je chí-kvadrát s J − ( p + 1) stupni volnosti. Pokud však J ≈ n , vzrůstá počet parametrů stejně rychle jako velikost souboru. Proto je p-value vypočtená pro tyto dvě statistiky v případě J ≈ n použitím rozdělení χ 2 ( J − p − 1) nesprávná. Jedním způsobem, jak se tomuto chybnému výpočtu vyhnout, je zavést určité seskupování. Na tom jsou založeny testy v následující kapitole.
3.5. Hosmer – Lemeshowy testy Hosmer a Lemeshow v [1] a [2] navrhli seskupování založené na hodnotách odhadnutých pravděpodobností. Předpokládejme J = n . Budeme předpokládat, že n sloupců odpovídá n hodnotám odhadnutých pravděpodobností, kde první sloupec náleží nejmenší hodnotě a n-tý sloupec největší hodnotě. Byly navrženy dvě seskupovací strategie. V první rozdělíme tabulku podle percentilů odhadnutých pravděpodobností a v druhé ji rozdělíme podle pevných hodnot odhadnutých pravděpodobností. V první metodě při použití deseti skupin bude první skupina obsahovat n1′ = n / 10 ′ = n /10 subjektů subjektů s nejmenší hodnotou odhadnuté pravděpodobnosti a poslední n10 s nejvyšší hodnotou odhadnuté pravděpodobnosti. U druhé metody s volbou desíti skupin dostaneme hraniční body s hodnotami k /10 , k = 1, 2,..., 9 , kde každá skupina obsahuje hodnoty odhadnutých pravděpodobností mezi dvěma hraničními body, tedy první skupina obsahuje subjekty, pro které je odhadnutá pravděpodobnost menší nebo rovna 0,1, druhá hodnoty mezi 0,1 a 0,2 a poslední obsahuje hodnoty vyšší než 0,9. Pro řádky s y = 1 dostaneme odhady očekávaných hodnot jako sumu odhadnutých pravděpodobností přes všechny subjekty skupiny. Pro řádky s y = 0 dostaneme odhad očekávaných hodnot součtem jedna minus odhadnutá pravděpodobnost pro všechny subjekty ve skupině. Ať už použijeme seskupování podle první či druhé metody, dostaneme statistiku Cˆ pro Hosmer – Lemeshowův test dobré shody výpočtem Pearsonovy chí-kvadrát statistiky z tabulky g × 2 pozorovaných a odhadnutých četností. Pro výpočet použijeme vzorec: 2
( o − nk′π k ) Cˆ = ∑ k , k =1 nk′ π k (1 − π k ) g
- 31 -
Verifikace modelu ___________________________________________________________________________ kde nk′ je celkový počet subjektů v k-té skupině, ck značí počet kovariačních schémat v k-tém ck
m jπˆ j . Rozdělení testové statistiky pro Cˆ je dobře ′ n j =1 k ck
decilu, ok = ∑ y j a π k = ∑ j =1
aproximováno chí-kvadrát rozdělením s g − 2 stupni volnosti, Cˆ ∼ χ 2 ( g − 2 ) . Podle [3] je doporučené používat seskupovací metodu založenou na percentilech odhadovaných pravděpodobností z důvodu větší podobnosti s χ 2 ( g − 2 ) . Rovněž se doporučuje používat deset skupin g = 10 . Tyto skupiny se pak nazývají decily rizik. V [4] autoři doporučují k výše uvedené metodě použít metodu normální aproximace distribuce Pearsonovy chí-kvadrát statistiky poprvé popsané Osiem a Rojekem v [5]. Postup je následující: 1. Získání a uložení vypočtených hodnot modelu πˆ j , j = 1, 2,..., J . 2. Vytvoření proměnné v j = m jπˆ j (1 − πˆ j ) , j = 1, 2,..., J . 3. Vytvoření proměnné c j =
1 − 2πˆ j vj
, j = 1, 2,..., J . 2
J
4. Vypočtení Pearsonovy chí-kvadrát statistiky Χ = ∑ j =1
(y
j
− m jπˆ j ) vj
2
.
5. Provedení vážené lineární regrese c na kovariátech modelu x s váhami v . Odtud dostaneme residuální součet čtverců a označíme jej RSS . J 1 6. Vypočtení opravy pro rozptyl A = 2 J − ∑ j =1 m j
.
Χ2 − J + p +1 . A + RSS 8. Vypočtení oboustranné p-value z normovaného normálního rozdělení.
7. Vypočtení statistiky z =
Výše uvedené testy nám poskytly metodiku, jak získat p-value. Tato obdržená hodnota nám umožní provést rozhodnutí pro takto stanovené hypotézy:
H0: Námi získaný model dobře popisuje naše data HA: Náš model popisuje data špatně a je třeba vytvořit jiný model přidáním či odebráním nezávislých proměnných, případně získáním většího množství dat
- 32 -
Popis algoritmů ___________________________________________________________________________
4. Popis algoritmů V této kapitole budou popsány vytvořené algoritmy, konkrétněji jejich požadovaný vstup a interpretace jejich výstupů.
4.1. R Všechny funkce byly napsány v programovacím jazyce a prostředí R. Jedná se o rozšíření komerčního produktu S, které je ovšem implementováno pod svobodnou licencí. Tento jazyk je hlavně využívaný pro statistickou analýzu a rovněž obsahuje funkce pro grafické zobrazení těchto dat. V základní podobě je R v podstatě příkazový řádek, existují ovšem i rozšíření přidávající grafické rozhraní (např. R Studio). Obliba tohota jazyka vede i jeho obsažení v jiných, někdy i komerčních, aplikacích, jakými jsou SPSS a Open Office. Popularita tohoto jazyka je především dána možností vytvářet balíčky obsahující funkce a procedury, které nejsou v základní verzi obsaženy. Tyto balíčky je možné stáhnout rovnou v základním rozhraní. Další výhodou R je možnost propojení s jazyky C, C++, Java, nebo Python. Všechny dále zmíněné funkce byly napsány ve verzi R 3.0.2 (2013-09-25) pro Windows x64.
4.2. Požadavky na vstup Všechny dále zmíněné funkce mají minimálně dvě vstupní hodnoty. Datový soubor a vektor parametrů logistické regrese β . R nativně nepodporuje čtení dat z excelovského souboru .xls, má ovšem zabudovanou funkci pro čtení souborů s příponou .csv a MS Excel převod do tohoto formátu přímo umožňuje. Pro načtení .csv souboru do R se pak používá funkce data=read.csv(„nazev_souboru.csv“,header=T,sep=“;“). Parametr header je binární a má hodnotu T, pokud otevíraný soubor obsahuje první řádek s popisky sloupců, jinak má hodnotu F. Parametr sep je pro oddělovač, který v .csv souboru odděluje jednotlivé sloupce. MS Excel při převodu z .xls na .csv používá jako oddělovač středník ;, při použití jiných metod je nutno - 33 -
Popis algoritmů ___________________________________________________________________________ brát na tuto položku zřetel. Samotný datový soubor pak musí obsahovat všechny nezávislé (vysvětlující) proměnné a jako poslední sloupec závislou (vysvětlovanou) proměnnou. Parametry logistické regrese β musí být na vstup přivedeny v podobě vektoru, nejjednodušeji příkazem beta=c( β1 ,…, β n , β 0 ). Pořadí jednotlivých koeficientů musí být stejné jako pořadí jednotlivých sloupců v datovém souboru, absolutní člen je uváděn jako poslední.
4.3. Lineární a exponenciální analýza Funkce, pro získání lineární a exponenciální analýzy jsou tři. Pro lineární analýzu je to linan.r, pro exponenciální analýzu expan.r a pro vylepšenou exponenciální analýzu enhexpan.r. Použitím příkazu source danou funkci načteme a pak ji můžeme zavolat. Hlavičky těchto funkcí: linan=function(betas,datas) expan=function(betas,datas) enhexpan=function(betas,datas) Všechny tři mají stejné vstupní parametry. První je vektor koeficientů β a druhým jsou načtená data. Výstupem je pak tabulka, která je vidět na obrázku obr. 1.
obr. 1: Výstup algoritmu pro exponenciální analýzu v R
- 34 -
Popis algoritmů ___________________________________________________________________________
4.4. Funkce vracející p-value Funkce, které vrací p-value, jsou celkem čtyři a to pearson.r používající Pearsonovu statistiku, resdev.r založená na výpočtu residuí deviance, osroj.r používající metodiku navrženou Osiem a Rojekem a hoslem.r podle postupu Hosmera a Lemeshowa. Hlavičky těchto funkcí: pearson=function(betas,datas,printcov=FALSE) resdev=function(betas,datas,printcov=FALSE) OsRoj=function(betas,datas,printcov=FALSE) HosLem=function(betas,datas,printcov=FALSE,groups=10) První vstupní parametr je opět pro hodnoty β , druhý pro vstupní data. Třetí parametr printcov volaný s hodnotou TRUE způsobí vypsání kovariačních schémat vstupních dat. Čtvrtý parametr u funkce HosLem slouží k volbě počtu skupin. Všechny tyto funkce vrací hodnotu p-value tak, jak je popsáno v teoretické části u jednotlivých metodik.
- 35 -
Zpracování dodaných dat ___________________________________________________________________________
5. Zpracování dodaných dat V této kapitole využijeme všech teoretických poznatků získaných v předchozích kapitolách. Použijeme příslušný software pro získání modelů pro dostupná data a využijeme vytvořených funkcí pro zhodnocení těchto modelů.
5.1. Datový soubor K otestování našich algoritmů nám byl poskytnut datový soubor z FNO. Tento soubor obsahuje údaje o operovaných pacientech s diagnózou rakoviny kolon a rekta. Rozsah dodaného souboru činil 666 pacientů, avšak čtyři záznamy musely být ze souboru vyřazeny kvůli chybějícím údajům. Tyto údaje byly sbírány postupně v letech 2001 – 2006. Základní charakteristiky souboru zachycují následující grafy.
Pohlaví
muž
290; 44% 372; 56%
žena
obr. 2: Koláčový graf pro rozložení pohlaví datového souboru
Z grafu je patrný vyšší počet pacientů mužského pohlaví. Tento nepoměr však pro nás neznamená žádné komplikace, skórovací systém pohlaví pacienta nezohledňuje a tato informace je zde uvedena pouze pro lepší ilustraci dodaného souboru.
- 36 -
Zpracování dodaných dat ___________________________________________________________________________
Metodika
304; 46%
otevřeně laparoskopicky
358; 54%
obr. 3: Koláčový graf pro způsob operace
Operační metodika je pro nás důležitý údaj. Model budeme vytvářet zvlášť pro pacienty operované laparoskopicky a zvlášť pro pacienty operované otevřenou metodou. Z grafu je patrné poměrně rovnoměrné rozložení obou skupin.
obr. 4: Krabicový graf zachycující vě ko vé rozložení pacientů
Na posledním grafu je zobrazeno věkové rozmezí pacientů. Nejmenší věk je 17 let, největší 97. Průměrný věk v souboru je 63,14 let, medián 64 let. Hodnota prvního resp. třetího kvartilu je 57, resp. 73 let a směrodatná odchylka má hodnotu 13 let.
- 37 -
Zpracování dodaných dat ___________________________________________________________________________
Mužů Žen Vyskytly Nevyskytly Minimum Maximum Průměr Směrodatná odchylka Medián První kvartil Třetí kvartil
Pohlaví Komplikace
V ěk
Laparoskopicky 214 144 151 207 18 97 62,3 13,2 64 57 71
Otevřeně 158 146 124 180 17 91 64,1 12,7 65 56,5 74
tab. 8: Rozložení veličin pro danou metodi ku
Protože budeme hledat modely odděleně pro laparoskopické a otevřené operace, je rozdělení základních veličin rozepsáno v tabulce.
5.1.1. Operační a fyziologické skóre, model POSSUM Již v úvodu byl zmíněn model POSSUM, který slouží pro odhad výskytu komplikací u pacientů s rakovinou kolorekta. Na počátku studie bylo 62 rizikových faktorů, které se postupně zredukovaly na 18 nezávislých faktorů. Tyto faktory se dále dělí na 12 faktorů popisujících fyziologický stav pacienta před operací a 6 faktorů samotného chirurgického výkonu. Tyto dílčí faktory jsou oceněny hodnotami 1,2,4 nebo 8, poté jsou tyto hodnoty sečteny a výsledek je označený jako fyziologické skóre pro 12 faktorů fyziologického stavu pacienta a operační skóre pro 6 faktorů operačního zákroku. Následující dvě tabulky převzaté z [6] a [7] popisují získání těchto hodnot. Operační skóre Závažnost a rozsah operačního výkonu Vícečetné operace (v posledních 30 dnech) Celková ztráta krve (ml) Kontaminace peritoneální dutiny Přítomnost malignity Naléhavost operace
Skóre 1
2
4
malý
střední
velký
1
2
>2 ≥ 1000
≤ 100
101-500
501-999
žádná
minimální (serózní)
lokálně hnis
žádná
jen primární
pozitivní uzliny naléhavá, je možná příprava > 2hod, operace do 24 hod od přijetí
elektivní
tab. 9: Zís kání operačního s kóre
- 38 -
8 komplexní, rozsáhlý
volný střevní obsah, hnis, krev vzdálené metastázy naléhavá, výkon je nutný do méně než 2 hod
Zpracování dodaných dat ___________________________________________________________________________ Skóre
Fyziologické skóre
1
2
4
Věk (roky)
≤ 60
≥ 71
Kardiální příznaky
bez selhávání
61-70 diuretika, digoxin, steroidy, terapie anginy pectoris nebo hypertenze
Rentgen srdce a plic Respirační příznaky
bez dušnosti
námahová dušnost
110-130
mírná CHOPCH 131-170 100-109 81-100 40-49
Rentgen plic Systolický krevní tlak (mm Hg) Tepová frekvence (minutová) Glasgow coma score Hemoglobin (g/l) 12
Leukocyty (.10 /l) Urea v séru (mmol/l) Natrium v séru (mmol/l) Kalium v séru (mmol/l)
Elektrokardiogram
periferní otoky, warfarin,
hraniční kardiomegalie hraniční dušnost (jedno patro) střední CHOPCH ≥ 171 90-99
8
zvýšený jugulární tlak
kardiomegalie klidová dušnost, (≥ 30/min) fibroza nebo konsolidace ≤ 89
101-120
≥ 121 ≤ 39
12-14
9-11
≤8
115-129 161-170 10,1-20,0 3,1-4,0
100-114 171-180 ≥ 20,1 ≤ 3,0
≤ 99 ≥ 181
≤ 7,5
7,6-10
10,1- 15,0
≥ 15,1
≥ 136
131-135
126-130
≤ 125
3,5-5,0
3,2-3,4 5,1-5,3
2,9-3,1 5,4-5,9
≤ 2,8 ≥ 6,0 jiný abnormální rytmus, ≥ 5 extrasystol /min, Q vlny nebo změny ST/T vlny
50-80 15 130-160 4-10
fibrilace síní (60-90/min)
normální
tab. 10: Zís ká ní f yziologic kého s kóre
Aplikací logistické regrese bylo riziko výskytu komplikace R vyjádřeno vztahem R ln = −5, 91 + 0,16 ⋅ FS + 0,19 ⋅ OS , kde FS značí fyziologické skóre a OS operační 1− R skóre.
5.2. Aplikace známých modelů Než vytvoříme model šitý na míru našim datům, zkusíme na tato data aplikovat již vytvořené modely. Prvním bude v předchozí kapitole zmíněný model POSSUM, druhým bude - 39 -
Zpracování dodaných dat ___________________________________________________________________________ model popsaný v [13]. Druhý model je ve tvaru R ln = −2, 75257 + 0, 08564 ⋅ FS + 0, 0654 ⋅ OS . 1− R Nejdříve aplikujeme algoritmy pro lineární, exponenciální a vylepšenou exponenciální analýzu popsané v kapitole 4.3. Lineární analýza POSSUM
Druhý model
Skupina [%]
Pacientů
Komplikací
Předpověď
Poměr
Pacientů
Komplikací
Předpověď
Poměr
<10 10-19 20-29 30-39 40-49 50-59 60-69 70-79 80-89 90-100 0-100
21 128 72 56 28 19 17 10 3 4 358
10 44 32 26 12 10 8 6 0 3 151
1 19 18 19 12 10 11 7 3 4 104
10 2,32 1,78 1,37 1 1 0,73 0,86 0 0,75 1,45
0 0 121 139 65 25 4 4 0 0 358
0 0 45 60 33 9 1 3 0 0 151
0 0 30 48 29 14 3 3 0 0 127
1 1 1,5 1,25 1,14 0,64 0,33 1 1 1 1,19
tab. 11: Lineární analýza pro pacienty operované laparoskopic ky
Lineární analýza POSSUM
Druhý model
Skupina [%]
Pacientů
Komplikací
Předpověď
Poměr
Pacientů
Komplikací
Předpověď
Poměr
<10 10-19 20-29 30-39 40-49 50-59 60-69 70-79 80-89 90-100 0-100
14 86 62 56 24 30 14 14 4 0 304
1 31 24 22 15 15 7 6 3 0 124
1 12 15 19 11 16 9 10 3 0 96
1 2,58 1,6 1,16 1,36 0,94 0,78 0,6 1 1 1,29
0 0 82 115 74 25 6 2 0 0 304
0 0 26 44 39 10 3 2 0 0 124
0 0 20 40 33 14 4 1 0 0 112
1 1 1,3 1,1 1,18 0,71 0,75 2 1 1 1,11
tab. 12: Lineární analýza pro pacienty operované otevřenou metodou
Tabulky lineární analýzy nám ukazují, že model POSSUM není moc vhodný pro předpověď komplikací. Pro oba soubory je předpovídaný počet komplikací výrazně nižší, než skutečný. Pro druhý model máme obdobnou situaci, ovšem rozdíl není natolik výrazný.
- 40 -
Zpracování dodaných dat ___________________________________________________________________________ Exponenciální analýza POSSUM
Druhý model
Skupina [%]
Pacientů
Komplikací
Předpověď
Poměr
Pacientů
Komplikací
Předpověď
Poměr
0-20 10-20 20-100 30-100 40-100 50-100 60-100 70-100 80-100 90-100 0-100
149 128 209 137 81 52 34 17 7 4 358
54 44 97 65 39 26 17 9 3 3 151
15 13 42 41 32 26 20 12 6 4 57
3,6 3,38 2,31 1,58 1,22 1 0,85 0,75 0,5 0,75 2,65
0 0 358 237 98 33 8 4 0 0 358
0 0 151 106 46 13 4 3 0 0 151
0 0 72 71 39 16 5 3 0 0 72
1 1 2,1 1,49 1,18 0,81 0,8 1 1 1 2,1
tab. 13: Exponenciální analýza pro pacienty operované laparos kopic ky
Exponenciální analýza POSSUM Skupina [%]
0-30 10-30 20-30 30-100 40-100 50-100 60-100 70-100 80-100 90-100
0-100
Druhý model
Pacientů
Komplikací
Předpověď
Poměr
162 148 62 142 86 61 32 18 4 0 304
56 55 24 68 46 30 16 9 3 0 124
24 15 12 43 34 30 19 13 3 0 67
2,33 3,67 2 1,58 1,35 1 0,84 0,69 1 1 1,85
Skupina [%]
0-10 10-20 20-30 30-100 40-100 50-100 60-100 70-100 80-100 90-100
0-100
Pacientů
Komplikací
Předpověď
Poměr
0 0 82 222 107 33 8 2 0 0 304
0 0 26 98 54 15 5 2 0 0 124
0 0 16 67 43 16 5 1 0 0 79
1 1 1,62 1,46 1,25 0,94 1 2 1 1 1,57
tab. 14: Exponenciální analýza pro pacienty operované otevřenou me todou
Situace u exponenciální analýzy je obdobná jako v případě lineární analýzy. V tomto případě je předpověď již výrazně horší, než-li je skutečnost. Ve srovnání modelu POSSUM a modelu, který byl vytvořen na geograficky obdobných datech, jaká máme k dispozici, vychází druhý model jako výhodnější. Tabulka pro vylepšenou exponenciální analýzu uvedena nebude, jelikož dvě předchozí metody hovoří výrazně proti modelům. Přece jen se ale jedná o slabší metody, proto ještě bude uvedena tabulka obsahující p-value, které získáme pomocí metod uvedených v kapitole 4.4.
- 41 -
Zpracování dodaných dat ___________________________________________________________________________
Residua deviancí Pearsonův test Hosmer-Lemeshow Osius-Rojek
POSSUM Laparoskopicky Otevřeně 0,00002334 0,0000005 0,00031172 0,00000379 0 0 0,02273075 0,00003586
Druhý model Laparoskopicky Otevřeně 0,00000518 0,00000167 0,1558699 0,07192259 0 0,00063562 0,12409096 0,00909371
tab. 15: P-value pro jednotlivé testy
V tab. 15 vidíme p-value pro jednotlivé typy testů. Zatímco model POSSUM můžeme zamítnout jako vhodný pro všechny provedené testy, druhý model se ukazuje jako poměrně vhodný pro použití na datech získaných z laparoskopických operací. Vidíme, že dva testy vykazují hodnoty blízké nule, ovšem další dva (Pearsonův a Osius-Rojek) mají přesvědčivé hodnoty pro použití modelu. Pro použití tohoto modelu hovoří i výsledek lineární analýzy, avšak exponenciální analýza neudává přesvědčivý výsledek. Rovněž počet kovariačních schémat u obou skupin operovaných pacientů je výrazně menší než celkový počet dat, měli bychom proto brát hodnoty prvních dvou testů jako významnější. Použití druhého modelu pro pacienty operované laparoskopicky tedy nemůžeme jednoznačně doporučit, ale ani zamítnout.
5.3. Vytvoření nového modelu Protože předchozí modely nepřinesly uspokojivý výsledek, pokusíme se o vytvoření nových modelů pro pacienty, kteří podstoupili zákrok laparoskopicky, a pro pacienty, kteří podstoupili otevřený zákrok. U obou těchto skupin známe fyziologické skóre a operační skóre, tyto dvě hodnoty budou sloužit jako nezávislé proměnné, a informaci o výskytu morbidity, která bude sloužit jako závislá (vysvětlovaná) proměnná. Hledaný model bude tedy R ve tvaru ln = β 0 + β1 ⋅ FS + β 2 ⋅ OS . 1− R
- 42 -
Zpracování dodaných dat ___________________________________________________________________________
obr. 5: Výs kyt ko mpli kací v závislosti na FS a OS u pacientů operovaných laparoskopic ky
obr. 6: Výs kyt ko mpli kací v závislosti na FS a OS u pacientů operovaných otevřenou me todou
Použitím statistického softwaru Statgraphics dostaneme modely: - 43 -
Zpracování dodaných dat ___________________________________________________________________________
R = −1, 0439 − 0, 00698 ⋅ FS + 0, 0628 ⋅ OS pro laparoskopický zákrok 1− R R ln = −2, 3812 + 0, 06869 ⋅ FS + 0, 06832 ⋅ OS pro otevřený zákrok. 1− R ln
Dále se podíváme na významnost jednotlivých koeficientů. Využijeme test věrohodnostního poměru, který je ve Statgraphicu implementován. Koeficient FS OS
χ2 0,0580942 5,52265
Stupně volnosti 1 1
p-value 0,8095 0,0188
Pro pacienty, kteří byli operováni laparoskopicky, vidíme, že na hladině významnosti 95% můžeme považovat koeficient u fyziologického skóre FS za nulový. Až budeme provádět verifikaci modelu, vyzkoušíme model, který tuto proměnnou obsahuje, a model, který ne.
Koeficient FS OS
χ2 5,467 5,15177
Stupně volnosti 1 1
p-value 0,0194 0,0232
U pacientů operovaných otevřenou metodou můžeme prohlásit, že na hladině významnosti 95% považujeme oba koeficienty jako významné.
5.4. Verifikace nového modelu Využijeme námi vytvořené algoritmy pro ověření vhodnosti nově vygenerovaných modelů. Pro laparoskopicky operované pacienty budeme uvažovat model obsahující veličinu FS a model bez této proměnné.
- 44 -
Zpracování dodaných dat ___________________________________________________________________________
5.4.1. Verifikace modelu pro laparoskopickou operaci Lineární analýza Model s FS
Model bez FS
Skupina [%]
Pacientů
Komplikací
Předpověď
Poměr
Pacientů
Komplikací
Předpověď
Poměr
<10 10-19 20-29 30-39 40-49 50-59 60-69 70-79 80-89 90-100 0-100
0 0 0 166 153 36 3 0 0 0 358
0 0 0 65 63 21 2 0 0 0 151
0 0 0 57 68 20 2 0 0 0 147
1 1 1 1,14 0,93 1,05 1 1 1 1 1,03
0 0 0 130 156 69 3 0 0 0 358
0 0 0 53 60 36 2 0 0 0 151
0 0 0 45 69 38 2 0 0 0 154
1 1 1 1,18 0,87 0,95 1 1 1 1 0,98
tab. 16: Lineární analýza pro pacienty operované laparoskopic ky
Podle lineární analýzy náš model predikuje velice dobře. Absence fyziologického skóre má za následek drobného zvýšení počtu predikovaných komplikací, což ale z praktického hlediska není na škodu. Exponenciální analýza Model s FS
Model bez FS
Skupina [%]
Pacientů
Komplikací
Předpověď
Poměr
Pacientů
Komplikací
Předpověď
Poměr
0-30 10-30 20-30 30-100 40-100 50-100 60-100 70-100 80-100 90-100 0-100
0 0 0 358 192 39 3 0 0 0 358
0 0 0 151 86 23 2 0 0 0 151
0 0 0 107 77 20 2 0 0 0 107
1 1 1 1,41 1,12 1,15 1 1 1 1 1,41
0 0 0 358 228 72 3 0 0 0 358
0 0 0 151 98 38 2 0 0 0 151
0 0 0 107 91 36 2 0 0 0 107
1 1 1 1,41 1,08 1,06 1 1 1 1 1,41
tab. 17: Exponenciální analýza pro pacienty operované laparos kopic ky
Exponenciální analýza vykazuje horší předpověď, komplikací určuje výrazně méně, než nastalo. Rozdíl mezi modelem obsahujícím fyziologické skóre a modelem bez něj je téměř nulový.
- 45 -
Zpracování dodaných dat ___________________________________________________________________________
Vylepšená exponenciální analýza Model s FS
Model bez FS
Skupina [%]
Pacientů
Komplikací
Předpověď
Poměr
Pacientů
Komplikací
Předpověď
Poměr
0-100 10-100 20-100 30-100 40-100 50-100 60-100 70-100 80-100 90-100 0-100
358 358 358 358 192 39 3 0 0 0 358
151 151 151 151 86 23 2 0 0 0 151
151 151 151 151 89 21 2 0 0 0 151
1 1 1 1 0,97 1,1 1 1 1 1 1
358 358 358 358 228 72 3 0 0 0 358
151 151 151 151 98 38 2 0 0 0 151
160 160 160 160 109 39 2 0 0 0 160
0,94 0,94 0,94 0,94 0,9 0,97 1 1 1 1 0,94
tab. 18: Vylepšená exponenciální analýza pro pacienty operované laparos kopic ky
Vylepšená exponenciální analýza ukazuje téměř přesnou shodu skutečných a predikovaných hodnot. Vynechání fyziologického skóre z modelu má za následek drobné zvýšení předpovídaných komplikací.
Residua deviancí Pearsonův test Hosmer-Lemeshow Osius-Rojek
Model s FS 0 0,0163 0 0,0000356
Model bez FS 0 0,00507 0 0,00000021
tab. 19: P-value pro jednotlivé testy
Ve skupině laparoskopicky operovaných pacientů máme celkem 358 záznamů, kovariačních schémat je 191. Pro tento model tedy neplatí J ≈ n a je vhodné používat spíše test na residua deviancí a Pearsonův test. Pokud se podíváme na test modelu, u kterého vyřadíme parametr fyziologického skóre, vidíme, že oba testy shodně zamítají vhodnost tohoto modelu. V případě, kdy fyziologické skóre zahrneme, test na residua model opět zamítá, ovšem hodnota Pearsonova testu na hladině významnosti 99% model nezamítá, pro hladinu významnosti 95% ale k zamítnutí již dochází. Nacházíme se tedy v nerozhodné oblasti. Výsledky lineární a exponenciální analýzy rovněž poukazují na vhodnost navrhovaného modelu. Protože náš model má za úkol pomoci lékařům při předpovídání operačních komplikací a v tomto případě očekáváme vyšší kvalitu modelu, je rozumné rozšířit datový soubor a vytvořit nový model na základě většího množství informací.
- 46 -
Zpracování dodaných dat ___________________________________________________________________________
5.4.2. Verifikace modelu pro otevřenou operaci Lineární analýza Skupina [%]
Pacientů
Komplikací
Předpověď
Poměr
<10 10-19 20-29 30-39 40-49 50-59 60-69 70-79 80-89 90-100 0-100
0 0 26 125 101 45 7 0 0 0 304
0 0 6 42 51 21 4 0 0 0 124
0 0 6 43 45 25 5 0 0 0 124
1 1 1 0,98 1,13 0,84 0,8 1 1 1 1
tab. 20: Lineární analýza pro pacienty operované otevřenou metodou
Tabulka ukazuje dobrou schopnost predikce našeho modelu. V jednotlivých procentuálních skupinách má tendenci nadhodnocovat riziko, to je ale ve výsledku pro pacienta příznivější. Exponenciální analýza Skupina [%]
Pacientů
Komplikací
Předpověď
Poměr
0-10 10-20 20-30 30-100 40-100 50-100 60-100 70-100 80-100 90-100 0-100
0 0 26 278 153 52 7 0 0 0 304
0 0 6 118 76 25 4 0 0 0 124
0 0 5 83 61 26 4 0 0 0 88
1 1 1,2 1,42 1,25 0,96 1 1 1 1 1,43
tab. 21: Exponenciální analýza pro pacienty operované otevřenou me todou
Použitá exponenciální analýza ukazuje poněkud horší schopnost predikce našeho modelu. Oproti skutečnosti má velmi optimistický náhled na zkoumaná data a to svědčí spíše proti testovanému modelu.
- 47 -
Zpracování dodaných dat ___________________________________________________________________________ Vylepšená exponenciální analýza Skupina [%]
Pacientů
Komplikací
Předpověď
Poměr
0-100 10-100 20-100 30-100 40-100 50-100 60-100 70-100 80-100 90-100 0-100
304 304 304 278 153 52 7 0 0 0 304
124 124 124 118 76 25 4 0 0 0 124
124 124 124 117 73 29 5 0 0 0 124
1 1 1 1,01 1,04 0,86 0,8 1 1 1 1
tab. 22: Vylepšená exponenciální analýza pro pacienty operované otevřenou metodou
Tato tabulka opět hovoří ve prospěch nového modelu. Ve většině procentuálních skupin dochází k velmi přesnému odhadu výskytu komplikací. Rovněž se projevuje snaha spíše pesimistického odhadu, což je opět výhoda pro pacienty. Na základě provedených jednoduchých analýz můžeme model spíše doporučit. Podstatné jsou ale testy, jejichž výsledky jsou prezentovány hned v následující tabulce. Residua deviancí Pearsonův test Hosmer-Lemeshow Osius-Rojek
0,00000035 0,06699 0 0,0034012
tab. 23: P-value pro jednotlivé testy
Pacientů operovaných otevřenou metodou jsme měli k dispozici 304 a pro tato data máme 196 kovariačních schémat. Stejně jako u laparoskopické techniky máme i zde výrazně méně kovariačních schémat než celkového počtu dat, proto je vhodné použití Pearsonova testu a testu na residua deviancí. Pearsonův test tento model nezamítá, i výsledky analýz mluví v jeho prospěch. Proti mluví hodnota testu residuí i Hosmer-Lemeshowův a OsiusRojekův test, použití posledních dvou zmíněných je ale vhodnější v situacích, kdy počet kovariačních schémat je blízký celkovému počtu dat.
- 48 -
Závěr ___________________________________________________________________________
6. Závěr Hlavním cílem této práce bylo vytvoření algoritmů pro posouzení vhodnosti modelů logistické regrese. Pro prvotní zhodnocení modelu byly vybrány metody lineární, exponenciální a vylepšené exponenciální analýzy. Dále byly vybrány a algoritmizovány testy verifikující daný model. Za tímto účelem byly vybrány čtyři testy a to Pearsonův test, test na residua deviancí, Hosmer-Lemeshův test a Osius-Rojekův test. První dva testy je vhodné použít v případě, kdy máme v datech výrazně nižší počet kovariačních schémat oproti původnímu počtu dat, druhé dva testy pak používáme, když je počet kovariačních schémat přibližně rovný počtu dat. Algoritmizace byla provedena v prostředí jazyka R. Tyto testy pak byly využity k otestování známých modelů na dodaných datech a z těchto dat byl vygenerovaný nový model, který byl rovněž otestován vytvořenými metodami. Došli jsme k závěru, že žádný z testovaných modelů neodpovídá požadovaným kvalitám a bylo doporučeno získání obsáhlejšího datového souboru a vytvoření nového modelu z těchto dat.
- 49 -
Literatura a reference ___________________________________________________________________________
7. Literatura a reference [1] Hosmer, D.W., and Lemeshow, S. A goodness-of-fit test for the multiple logistic regression model. Communications in Statistics, 1980, A10, 1043-1069. [2] Lemeshow, S., and Hosmer, D.W. The use of goodness-of-fit statistics in the development of logistic regression models. American Journal of Epidemiology, 1982, 115, 92-106. [3] Hosmer, D.W., Lemeshow, S., and Klar, J. Goodness-of-fit testing for multiple logistic regression analysis when the estimated probabilities are small. Biometrical Journal, 1988, 30, 911-924. [4] Hosmer, D.W., Hosmer, T., Le Cessie, S., and Lemeshow, S. A comparison of goodnessof-fit tests for the logistic regression model, Statistics in medicine, 1997, 16, 965-980. [5] Osius, G., and Rojek, D. Normal goodness-of-fit tests for multinomial models with large degrees-of-freedom. Journal of the American statistical association, 1992, 87, 1145-1152. [6] Martínek, L. Aplikace specializovaných skórovacích systémů pro objektivizaci rizik laparoskopických operací kolorekta. Ostrava, 2006. Doktorská disertační práce. Chirurgická klinika Fakultní nemocnice s poliklinikou Ostrava. [7] Rabasová, M. Diskriminační analýza jako nástroj pro hodnocení chirurgických rizik. Ostrava, 2012. Doktorská dizertační práce. VŠB - Technická univerzita Ostrava [8] Hosmer, D.W., Lemeshow, S. Applied logistic regression, Wiley 2000, ISBN 0-47135632-8. [9] Wijesinghe, L.D., Mahmood, T., Scott, D.J.A., Berridge, D.C., Kent, P.J., Kester, R.C. Comparison of POSSUM and the Portsmouth predictor equation for predicting death following vascular surgery. Br. J. Surg., 1998, 85, s.209-212 [10] Jahoda P., Briš R., Korekce exponenciální analýzy při posuzování vhodnosti modelu morbidity, Sborník z konference REQUEST 2007, CQR 2007, ISBN 978-80-01-03709-6, str. 124-135. [11] R.Briš, M.Litschmannová, Statistika I, elektronické skriptum VŠB TUO, FEI,2004. [12] Briš R., Litschmannová M., STATISTIKA II., E-learningový prvek pro podporu výuky odborných a technických předmětů,v rámci projektu CZ.O4.01.3/3.2.15.2/0326, VŠB TU Ostrava, 2007, ISBN 978-80-248-1482-7.
- 50 -
Literatura a reference ___________________________________________________________________________ [13] Briš R.: Modelování rizika morbidity laparoskopických operací pomocí logistické regrese, 8. národní konference Statistické dny v Brně, 27. - 28. června, 2006, ISBN 80-214-3214-4. Pg.13-14.
- 51 -