Letný semester 2024 | Cvičenie 2 | 04.03.2024
Prihlásenie k SAS OnDemand: https://www.sas.com/en_us/software/on-demand-for-academics.html
Doporučená literatúra a ďalšie užitočné materiályPodobne ako v prípade štandardných (lineárnych) regresných modelov (resp. akejkoľvek štatistickej analýzy – t.j. aplikácie nejakého konkrétneho stochastického modelu), je štatistická analýza dat rozdelená do dvoch samostatných, ale vzájomne prepojených časti – exploratívnej analýzy a tzv. “confirmatory” analýzy. Exploratívna analýza využíva rôzne popisné charakteristiky (t.j., empirické odhady teoretických charakteristík) a nástroje pre vizualizáciu dat. Konfirmačná analýza (i.e., štatistická inferencia) hľadá pomocou vhodných pravdepodobnostných metód odpoveď na konkrétnu (expertnú) otázku (t.j., testuje konkrétnu nulovú hypotézu).
Hlavnou výhodou longitudinálnych dat (a modelov určených na ich analýzu) je to, že na rozdiel od cross-sekčného porovnania dvoch subpopulácii umožňujú odhadovať a porovnávať (testovať) aj zmenu v rámci konkrétnych subjektov a to vzhľadom k časovo sa meniacim sa hodnotám nezávislých premenných v rámci daného subjektu. Klasická regresná analýza cross-sekčných dat umožňuje len porovnávať dve podskupiny, ktoré sa vzájomne líšia vzhľadom k hodnotám nezávislých premenných. Longitudinálne data a modely teda umožňujuú porovnanie jednak statických (medzi-populačných) rozdielov a tiež dynamických zmien v rámci konkrétneho subjektu. Nejedná sa ale o kauzálnu analýzu.
Pre \(N \in \mathbb{N}\) nezavislých
subjektov, z ktorých každý subjekt \(i \in
\{1, \dots, N\}\) sledujeme opakovane \(n_i \in \mathbb{N}\) krát, môžeme data
reprezentovať v tvare \(\{(Y_{ij},
\boldsymbol{X}_{ij}^\top)^\top;~i = 1, \dots, N; j = 1, \dots,
n_{i}\}\). V prípade reprezentácie v zmysle \[
\{(\boldsymbol{Y}_{i}^\top, \boldsymbol{X}_{i 1}^\top, \dots,
\boldsymbol{X}_{i n_i}^\top)^\top;~i = 1, \dots, N\}
\] môžeme dokonca (v určitom zmysle) hovoriť o náhodnom výbere
(nezávislé a rovnako (združene) rozdelené náhodné vektory). Je nutn0 si
ale uvedomiť, že rozmer (dĺžka) náhodných vektorov \((\boldsymbol{Y}_{i}^\top, \boldsymbol{X}_{i
1}^\top, \dots, \boldsymbol{X}_{i n_i}^\top)^\top\) je pre rôzne
\(i \in \{1, \dots, N\}\) rôzna.
Obyčajný linárny regresný (tzv. “cross-sectional”) model môžeme formulovať ako \[ Y_{i 1} = \boldsymbol{x}_{i 1}^\top \boldsymbol{\beta}_{CS} + \varepsilon_{i 1}, \qquad i = 1, \dots, N, \] pre vektor neznámych parametrov \(\boldsymbol{\beta}_{CS} \in \mathbb{R}^p\). Príslušný ohad \(\widehat{\boldsymbol{\beta}}_{CS}\) lze po jednotlivých zložkách interpretovať ako odhadovaný očakávaný rozdiel medzi dvoma podskupinami, ak je rozdiel príslušného regresoru (jeden konkrétny element v \(\boldsymbol{x}_{i 1}\)) jednotkový a ostatné hodnoty regresorov (zvyšných elementov \(\boldsymbol{x}_{i 1}\)) sú rovnaké.
Na druhej strane, jednoduchý model pre longitudinálne data môžeme zapísať aj v tvare \[ Y_{i j} = \boldsymbol{x}_{i 1}^\top \boldsymbol{\beta}_{CS} + (\boldsymbol{x}_{i j} - \boldsymbol{x}_{i 1})^\top \boldsymbol{\beta}_{L} + \varepsilon_{i j}. \qquad i = 1, \dots, N; j = 1, \dots, n_i, \] kde celkový počet pozorovaní je \(\sum_{i = 1}^N n_{i}\). Model ale okrem odhadu pre vektorový parameter \(\boldsymbol{\beta}%/CS(\) (ktorý má ekvivalentnú interpretáciu ako v prípade obýčajného regresného modelu vyššie) umožňuje odhadnúť aj vektorový parameter \(\boldsymbol{\beta}_L\). Interpretácia tohto parametru je ale už specifická pre daný subjekt a vyjadruje zmenu v rámci daného subjektu \(i \in \{1, \dots, N\}\). Príslušný odhad \(\widehat{\boldsymbol{\beta}}_{L}\) môžeme po jednotlivých zložkách interpretovať ako očakávanú zmenu v rámci daného subjektu v čase, ak sa príslušný regresor zmení v čase o jednotku.
Pre ilustráciu použijeme datový súbor
pig_weights.csv, ktorý udáva váhu u 48
ošípaných, ktoré boli postupne merané po dobu 9 týždňov
(t.j.,
follow-up obdobie, po ktoré subjekty sledujeme, je 9 týždňov – nie ale
nutné stejných 9 týždňov). Do programu SAS data načítame následujúcim
príkazom:
libname weight '/home/u63241636/sasuser.v94';
filename reffile '/home/u63241636/sasuser.v94/data/pig_weights.csv';
proc import datafile=reffile
dbms=csv
out=weight.data
replace;
getnames=yes;
run;
proc print datafile = weight.data;
run;Samostatne sa podívajte na niektoré základné popisné charakteristiky spočítané z dat. Pre ilustráciu niektorých popisných charakteristík poslúži napr. aj následujúci graf. Pokúste sa ale vytvoriť aj niekoľko ďalších grafov a obrázkov samostatne.
proc sgscatter data=weight.data;
matrix weight week /
ellipse=(type=mean)
diagonal=(histogram kernel);
run;V nasledujúcej časti sa špecificky zameriame na vizualizáciu datového
súboru vhľadom k jeho časovo-zavislej štruktúre a tzv.
subject-specific longitudinálnych profilov.
Porovnajte následujúce grafické výstupy:
proc sgscatter data=weight.data;
plot weight * week ;
run;proc sgscatter data=weight.data;
plot weight * week / group = id;
run;proc sgplot data=weight.data;
series y=weight x=week /group=id;
run;Podívajte sa na základnú syntax procedrúry
PROC sgscatter napr. na
tejto
stránke a analogicky tiež na syntax procedúry
PROC sgplot napr. na
tejto
stránke.
Jednoduché longitudinálne profily nie vždy dobre poslúžia k tomu, aby
sme získali dôležité informácie o (časovom) vývoji konkrétneho
jednotlivca vzhľadom k celkovej skupine. Špeciálne v prípadoch, keď je
uvažovných subjektov príliš mnoho, je niekedy užitočné zvážiť aj graf s
tzv. štandardizovanými longitudinálnymi profilmi.
Jedna z možnosti, ako štandardizované longitudinálne profily získať v programe SAS, je v niekoľkých krokoch uvedená v následujúcich bodoch:
K vytvoreniu štandardizovaných profilov potrebujeme v prvom rade vhodné popisné charakteristiky spočítané z dat. Konkrétne výberové priemery a smerodatné odchýlky v jednotlivých (časových) meraniach. Vytvorenie popisných charakteristík (agregované data vzhľadom k jednotlivým týždňom follow-up obdobia) získame napr. následujúcim spôsobom:
proc means data=weight.data mean std;
class week;
var weight;
output out=summary mean = MeanValue std = StdValue;
run;proc print datafile = summary;
run;V ďalšom kroku potrebujeme priradit jednotlivé agregované (empirické) charakteristiky k jednotlivým pozorovaniam (spojiť dva datove súbory). V programe SAS je pre tento krok kľôčove správne usporiadanie jednolivých pozorovaní – v tomto konkrétnom prípade ide o usporiadanie oboch datových súborov podľa premennej, ktorá identifikuje jednotlivé týždne (premenná week):
proc sort data=weight.data; by week; run;
proc sort data=summary; by week; run;Zlúčenie oboch datových súborov dohromady a opätovne preusporiadanie podľa jednotlivých subjektov (ako v pôvodných datach):
data weight.dataStd;
merge summary weight.data ;
by week;
run;
proc sort data=weight.dataStd; by id; run;V aktuálnom datovom súbore už celkmo jednoducho (pomocou lineárnej kombinácie existujúcich premenných) vytvoríme novú premennú, ktorá bude obsahovať štandardizované profily: :
data weight.dataStd;
set weight.dataStd;
weightStd = (weight - MeanValue)/StdValue;
run;proc print datafile = weight.dataStd;
run;Momentálne sú data pripravené k vykresleniu štandardizovaných profilov. Jednotlivé “subject-specific profily” sú často výraznejšie práve v takomto type grafov:
proc sgplot data=weight.dataStd;;
series y=weightStd x=week /group=id;
run;Graf štandardizovaných rezíduí je možné napríklad použíť k nahliadnutiu na rezídua, prípadne na overenie niektorých predpokladov súvisiacich s náhodnými chybami.
ods graphics on;
proc loess data=weight.dataStd;
model weightStd=week;
run;V tejto časti použijeme trochu komplexnejší datový súbor o pacientoch so sklerózou multiplex (roztroušená skleróza). Podkladový datový súbor je možné stiahnúť ako csv súbor (následne je potrebné súbor uploadovať do SAS OnDemand a pomocou nasledujúceho SAS kódu načítať).
libname sm '/home/u63241636/sasuser.v94';
filename reffile '/home/u63241636/sasuser.v94/data/sm_data2.csv';
proc import datafile=reffile
dbms=csv
out=sm.data
replace;
getnames=yes;
run;
proc print datafile = sm.data;
run;V praxi sa niektoré premenné rozlíšujú do dvoch skupín: tzv. subject-specific variables a tzv. time-dependent variables. V prvom prípade sa jedná o premenné, ktoré sú v čase nemenné (napr. pohlavie pacienta, vek pacienta v čase podstúpenia liečby, zdravotný stav pacienta v čase podstúpenia liečby, alebo aj napr. výskyt/nevýskyt komplikácii počas sledovania). V tom druhom prípade sa jedná o premenné, ktoré sa vrámci konkrétneho jedinca v priebehu opakovaných meraní môžu meniť (napr. hodnota EDSS, alebo opäť výskyt komplikácii zaznamenaný samostatne pri každej kontrole).
Ide o zisťovanie a vyšetrovanie vzťahov (závislostnej štruktúry) v tzv. subject-specific premenných – porovnávanie vrámci nezávislých subjektov (podskupín, sub-populácii). V prípade datového súboru s pacientami so sklerózou multiplex – napríklad EDSS (t.j., Expanded Disability Status Scale) v závislosti na pohlaví pacienta:
proc loess data = sm.data;
model EDSS = age / clm smooth=0.1 0.5 residual;;
by gender;
run;Jedná sa o zisťovanie vplyvu/zmeny vrámci opakovaných pozorovaní – korelované pozorovania vrámci daného subjektu. Opäť v prípade datového súboru s pacientami so sklerózou multiplex – napríklad EDSS (a jeho vývoj) v závislosti na čase:
proc loess data = sm.data;
model EDSS = time / clm smooth=0.1 0.5 residual;
run;Nemusí sa ale nutne jednať pouze o čas ako takový – t.j., premennú time. Analogicky móžeme zisťovať napr. zmenu v EDSS v závislosti na zmene hodnoty v premennej relapse:
proc sort data=sm.data; by relapse; run;
proc boxplot data=sm.data;
plot EDSS*relapse;
run;Na záver sa opäť vrátime k datovému súboru o váhach ošípaných, ktorý
sme navyše doplnili o štandardizované profily. Vytvorené štandardizované
longitudinálne profily použijeme k tomu, aby sme nahliadli do
časovo-závislej štruktúry.
proc print datafile = weight.dataStd;
run;Existuje samozrejme množstvo rôznych spôsobov, ako empiricky
vyšetrovať a zisťovať prípadnú závislosť medzi pozorovaniami.
Predpokládame, že jednotlivé súbjekty sú vzájomne nezávisle (malo by to
vyplývať z povahy samotného experimentu). Opakované pozorovania vrámci
konkretného jedinca ale nezávislé samozrejme nie sú. Pre následné
modelovanie (t.j., “confirmatory analysis”) je ale dôležité závislostnú
štruktúru opakovaných pozorovaní vyšetriť.
V prvom kroku vhodne pripravíme data. Pri analýze longitudinálnych
dat sa často používa terminológia, ktorá odkazuje na tzv.
long-data štruktúru dat, prípadne tzv. wide-data typ dat.
Data uložené v datovom súbore weight.dataStd sú reprezentované v
tzv. long-data type. Pomocou PROC transpose ale data
môžeme reprezentovať vo wide-data formáte:
proc sort data=weight.dataStd; by id; run;
proc transpose data=weight.dataStd out=weight.dataWide;
by id;
var weightStd;
run;Wide-data formát už môžeme priamo využiť pre spočítanie variančenej kovariačnej (resp. korelačnej) matice
proc corr data=weight.dataWide;
var COL1 COL2 COL3 COL4 COL5 COL6 COL7 COL8 COL9;
run;prípadne tzv. matice scatterplotov, ktorá časovo-závislostnú štruktúru reprezentuje vizuálne:
proc sgscatter data=weight.dataWide;
title "Scatterplot Matrix -- Repeated observations";
matrix COL1 COL2 COL3 COL4 COL5 COL6 COL7 COL8 COL9;
run;
title;Predchádzajúci výstup lze také porovnať s následujúcim výstupom:
proc sgscatter data=weight.dataWide;
matrix COL1 COL2 COL3 COL4 COL5 COL6 COL7 COL8 COL9 / diagonal=(histogram kernel);
run;Pripomeňme, že pre ľubovolnú náhodnú veličinu \(X \sim F\) z nejakého (absolútne) spojitého rozdelenia s distribučnou funkciou \(F\) platí, že \[ P[X \leq x] = F(x), \qquad \forall x \in \mathbb{R}. \] Zároveň ale dostávame, že platí tiež \[ P[F(X) \leq x] = P[X \leq F^{-1}(x)] = F\Big(F^{-1}(x)\Big) = x, \qquad \textrm{pre} \quad x \in [0,1], \] kde \(F^{-1}(u)\) pre \(u \in [0,1]\) je príslučná kvantilová funkcia náhodnej veličiny \(X\) (t.j., inverzná funkcia k distribučnej funkcii \(F\)). Inými slovami, transformovaná náhodná veličina \(F(X)\) má rovnomerné (spojité) rozdelenie na intervale \([0,1]\).