6.6 Alkalmazások
1) 6.1. Táblázat nyújt az eredmények egy kohorsz tanulmány célja, hogy meghatározza, többek között, a használatának előnyei, mint egy szűrővizsgálat, terheléses teszt (EST) intézkedések. Ilyen vizsgálatok során a koszorúér-betegségek (CAD) diagnosztizálása során pass/fail eredmény származhat .
6.1. táblázat. Screening test for physical stress
| CAD | ||||
|---|---|---|---|---|
| Healthy | Ill | Total | ||
| EST | Negative | 327 | 208 | 535 |
| Positive | 115 | 815 | 930 | |
| Total | 442 | 1 023 | 1 465 | |
feltételezzük, hogy nincs ellenőrzési torzítás:
A)
ezen zavartmátrix alapján jelölje meg a következő értékeket (95% – os konfidencia intervallumokkal): érzékenység és specificitás, pozitív és negatív prediktív érték.
b)
mi a görbe alatti terület értéke az aktuális adatokhoz?
az adattábla létrehozása az alábbiak szerint történhet, a() :
> tab < – as.táblázat(mátrix(c(815,115,208,327), nrow = 2, byrow = IGAZ,
dimnames = lista(EST = c(“+”,”−”),
CAD = c(“+”,”-“))))
> lapon
Meg kell jegyezni, hogy a táblázat átszervezték, így megfelelnek a szokásos jelöléseket, a “pozitív” események az első sorban (jellemzően az expozíció), valamint az első oszlop (általában a betegség) a táblázat.
ezután az összes szükséges mennyiséget nagyon kevés aritmetikai művelettel lehet kiszámítani R-vel. Az epiR csomag azonban tartalmazza a járványügyi kérdésekre adott válaszhoz szükséges összes eszközt a 2 × 2 táblázatban. Így a parancs epi.a tesztek() prevalenciát, érzékenységet/specificitást, valamint pozitív és negatív prediktív értékeket adnak:
> könyvtár(epiR)
> epi.tesztek (tab)
a görbe alatti területet illetően külső parancs is alkalmazható, roc.- tól.táblázat(), a csomagban epicalc:
> könyvtár (epicalc)
> roc.- tól.táblázat (tabulátor, grafikon = hamis)
2) a 6.2 .táblázat összefoglalja a 40-59 éves férfiak körében megfigyelt myocardialis infarktusok arányát, akik esetében a vérnyomás szintjét és a koleszterinszintet mérték, rendezett osztályok formájában.
6.2. táblázat. Infarktus és vérnyomás
| koleszterin (mg/100 ml) | ||||||||
|---|---|---|---|---|---|---|---|---|
| BP | & / div > lt; 200 | 200 – 209 | > 284 | |||||
| & lt; 117 | 2/53 | 0/21 | 0/15 | 0/20 | 0/14 | 1/22 | 0/11 | |
| 117 – 126 | 0/66 | 2/27 | 1/25 | 8/69 | 0/24 | 5/22 | 1/19 | |
| 127 – 136 | 2/59 | 0/34 | 2/21 | 2/83 | 0/33 | 2/26 | 4/28 | |
| 137 – 146 | 1/65 | 0/19 | 0/26 | 6/81 | 3/23 | 2/34 | 4/23 | |
| 147 – 156 | 2/37 | 0/16 | 0/6 | 3/29 | 2/19 | 4/16 | 1/16 | |
| 157 – 166 | 1/13 | 0/10 | 0/11 | 1/15 | 0/11 | 2/13 | 4/12 | |
| 167 – 186 | 3/21 | 0/5 | 0/11 | 2/27 | 2/5 | 6/16 | 3/14 | |
| > 186 | 1/5 | 0/1 | 3/6 | 1/10 | 1/7 | 1/7 | 1/7 | |
The data are available in the file hdis.dat formájában egy táblázat négy oszlop jelzi, illetve a vérnyomás (nyolc kategória, eddig 1-8), az arány a koleszterin (hét kategória, eddig 1-7), a szám a miokardiális infarktus, valamint az összes egyének száma. Az egyesület között a vérnyomás, valamint a valószínűsége, hogy a szenvedés egy miokardiális infarktus az érdekes ponthoz:
a)
ki kell Számítani az arányokat miokardiális infarktus minden szint, a vérnyomás, valamint képviseli őket egy asztalnál grafikus formában;
b)
Express az arányok kiszámítása (a) a logit formában;
c)
Alapuló logisztikai regressziós modell határozza meg, hogy van-e szignifikáns egyesület a küszöbön α = 0.05 között a vérnyomás, kezelni, mint egy mennyiségi változót figyelembe véve az osztály központok, valamint a valószínűsége, hogy a szívroham;
d)
Kijelző ugyanazon a grafikon az empirikus arányok, valamint a logisztikai regressziós görbe szerint a vérnyomás értékek (osztály központok).
Az adatok importálása hdis.a dat a következőképpen érhető el:
> bp < – read.táblázat (“hdis.dat”, header = TRUE)
> str(bp)
mivel ellenőrizhető, nincs címke a kamatváltozók szintjéhez (bpress a vérnyomáshoz, chol a koleszterinszinthez). A címkék létrehozásához és társításához a következő parancsokat lehet beilleszteni:
> blab < – c(“<117”,”117-126”,”127-136”,”137-146”,
“147-156”,”157-166”,”167-186”,”>186”)
> clab < – c(“<200”,”200-209”,”210-219”,”220-244”,
“245-259”,”260-284”,”>284”)
> bp$bpress < – factor(bp$bpress, labels = blab)
> bp$chol < – factor(bp$chol, labels = clab)
az utolsó parancs az eredeti változókat kvalitatív változókká és társaikká alakítja át a blab és clab által definiált címkékre. Annak ellenőrzéséhez, hogy az adatbázis most a kívánt formában van-e, a következő parancsok hasznosak:
> str (bp)
> összefoglaló (bp)
a nyilatkozatban kitett relatív frekvenciák táblázata most reprodukálható:
Mivel mi fog összpontosítani a kapcsolat vérnyomás, szívinfarktus csak, szükséges, hogy az összesített adatokat a koleszterin szintet. Más szavakkal, az egyes vérnyomásszintek értékeit fel kell felhalmozni, függetlenül a koleszterinszinttől. Lehetőséget kell adni arra is,hogy átnevezzük a bpress változó szintjeit az osztályközök központjaival:
> blab2 < – c(111.5,121.5, 131.5,141.5,151.5,161.5,176.5,191.5)
> bp$bpress < – rep(blab2, minden = 7)
> dfrm < – összesített(bp,
lista(bpress = bp), sum)
Ez az utolsó parancs, összesített(), amely lehetővé teszi, hogy összesített adatok: az összegeket minden értékek (összeg) a változó choi, amelyek nem szerepelnek a listán változók őrizni elemzés. Ugyanakkor az eredményeket egy új dfrm nevű adatbázisban tárolják. Ne feledje, hogy ebben az esetben nem a képlet jelölését használjuk, hanem az aggregátum () alternatív szintaxisát (válaszváltozó, amelyet egy vagy több magyarázó változóból álló lista követ). Az összesített adatok áttekintése a head() bemenetként történő megadásával érhető el:
> head(dfrm, 5)
ha rendelkezésre áll egy P arány, annak értéke olyan skálán, amelynek egységeit a relációs napló adja(p/(1 − p))). A következő parancsokból származtatható az arányok konvertálása, úgy számítva, mint a hdis / total a logit egységekben:
> dfrm$logit < – logit(dfrm$hdi-k/dfrm$összesen)
megjegyzendő, hogy egy kis funkció került meghatározásra, mely értékek x, ami itt feltételezzük, hogy egy arányok, az annak megfelelő log(x/(1 − x)). Ugyanígy írhatnánk:
> log (((dfrm$hdis/dfrm$total)/(1-dfrm$hdis/dfrm$total))
a számított értékek közvetlenül a prop és logit név alatt kerültek az adatkeretbe.
A logisztikus regressziós modell van írva, a következő módon:
> glm(cbind(hdi-k,a teljes-hdi-k) bpress,data = dfrm,családi = binomiális)
A készítményben használják, cbind(hdi-k, a teljes-hdi-k) bpress, figyelembe veszi azt a tényt, hogy mi elérése csoportosított adatokat, nem pedig az egyes válaszok. A GLM() parancs a family = binomial opcióval egy logisztikai regressziónak felel meg, amely alapértelmezés szerint a Logit funkciót kanonikus linkként használja.
a regressziós együtthatók az összefoglaló () használatával jeleníthetők meg:
> összefoglaló (glm (cbind (hdis, total-hdis) bpress, data = dfrm,
család = binomiális))
az előző eredmény tartalmazza az alapvető információkat a vérnyomás és a szívinfarktus valószínűsége közötti összefüggés statisztikai jelentőségének megválaszolásához: a regressziós együttható (a log-odds skálán) 0,024, a szokásos 5% – os küszöbértéknél pedig jelentős (lásd a PR oszlopot (> / z/)).
a vizsgált vérnyomás különböző szintjeinek megfelelő szívroham valószínűsége a következő:
> glm.res < – glm(cbind(hdis, total-hdis) bpress, data = dfrm,
family = binomial)
> predict(glm.res)
meg kell jegyezni, hogy az R által generált közbenső eredményeket glm névvel tárolták.res a predict () parancs használata előtt. Az R által generált előrejelzéseket logit formában fejezzük ki, és összehasonlíthatjuk a megfigyelt és előre jelzett naplókat:
> cbind(dfrm, logit.predit = predict(glm.res))
gyakorlatilag minden elem rendelkezésre áll, hogy képviselje a megfigyelt és grafikusan előre jelzett myocardialis infarctus arányát a vérnyomás szintje szerint. A modell előrejelzései arányokban kifejezve, nem pedig a logodds skálán hiányoznak. Ezért kívánatos lehet az előrejelzési görbe rajzolása, vagyis a vérnyomás szerinti szívroham valószínűsége, anélkül, hogy ez utóbbit a bpress változó nyolc értékére korlátoznánk. Itt következik egy lehetséges megoldás:
> dfrm$prop.predit < – predict(glm.res, type = “response”)
> f < – function(x)
1/(1 + exp(-(coef(glm.res) + coef(glm.res)*x)))
> xyplot(hdis/total ˜ bpress, data = dfrm, aspect = 1.2, cex = .8,
xlab = “Blood pressure”,ylab = “Infarction likelihood”,
panel = function(x, y, …) {
panel.xyplot(x, y, col = “gray30”, pch = 19, …)
panel.curve(f, lty = 3, col = “gray70”)
panel.points(x, dfrm$prop.predit, col = “gray70”,…)
})
3) egy eset-ellenőrzési vizsgálat az alkoholfogyasztás és a dohány és a nyelőcsőrák közötti összefüggésre összpontosított emberekben (“Ille és Villaine” tanulmány). Az esetek csoportjába 200, nyelőcsőrákban szenvedő beteg tartozott, akiket 1972 januárja és 1974 áprilisa között diagnosztizáltak. Összesen 775 ellenőrző férfit választottak ki a választási listáról. 6. táblázat.3 mutatja az alanyok összességének eloszlását napi alkoholfogyasztásuk szerint, figyelembe véve, hogy a 80 g-nál nagyobb fogyasztást kockázati tényezőnek tekintik .
6.3. táblázat. Alkoholfogyasztás és nyelőcsőrák
| alkoholfogyasztás (g/nap) | |||
|---|---|---|---|
| ≥ 80 | Teljes | ||
| Case | 96 | 104 | 200 |
| Control | 109 | 666 | 775 |
| Teljes | 205 | 770 | 975 |
a)
Mi az értéke annak az esélye, az arány pedig a 95% – os konfidencia-intervallum (Woolf módszer)? Jó becslés a relatív kockázatról?
b)
a veszélyeztetett fogyasztók aránya azonos az esetek és a kontrollok között (lásd α = 0, 05)?
c)
készítse el a logisztikus regressziós modellt, amely lehetővé teszi az alkoholfogyasztás és az egyének státusza közötti összefüggés tesztelését. A regressziós együttható jelentős?
d)
keresse meg a megfigyelt esélyarány B) pont szerinti értékét és annak konfidenciaintervallumát a regressziós analízis eredményei alapján.
mivel az (egyedi) nyers adatok nem állnak rendelkezésre, közvetlenül a problémamegállapításban megadott frekvenciatáblázattal kell dolgozni:
> alcohol < – matrix(c(666,104,109,96), nr = 2,
dimnames = list(c(“Control”,”Case”),
c(“<> = 80”)))
> alcohol
Regarding the odds ratio, the command oddsratio() from the package vcd will be used:
> library(vcd)
> oddsratio(alkohol, log = FALSE)
a log = FALSE opció biztosítja, hogy a visszaadott eredmény helyesen megfeleljen az esély aránynak, nem pedig az esély Arány naplójának. Az aszimptotikus konfidencia intervallum eléréséhez a következőket kell alkalmazni:
> confint(oddsratio(alkohol, log = hamis))
hasonlóképpen összefoglalót(oddsratio(alkohol)) is lehet használni a log odds Arány hipotézisének elvégzésére (H0 : log(θ) = 0).
a parancs prop.teszt () használható annak a hipotézisnek a tesztelésére, hogy a napi fogyasztással rendelkező személyek aránya ≥ 80 g, és azonos az esetekben, valamint a kontrollokban, jelezve a problémamegállapításban megadott kereszttáblázatból megfigyelt értékeket:
> prop.teszt (c (96,109), c(200,775), correct = FALSE)
Ez a teszt pontosan megegyezik egy Z-teszttel, amely az adatokból becsült két arány közötti különbséget teszteli (ha nem alkalmaznak folytonossági korrekciót).
a logisztikai regressziós modell tekintetében a készenléti táblázatot olyan adattáblává kell alakítani, amelyben a két minőségi változó (betegség és expozíció) megfelelően van ábrázolva, azaz adatkerettel. Ezt a reshape2 csomagból származó melt() paranccsal lehet elvégezni, amely lehetővé teszi az adatok táblázatos formában “hosszú” formátumba történő átalakítását. Ebben a szakaszban a betegség változó szintjét 0 (kontroll) és 1 (eset) értékkel kell újradefiniálni, bár ez nem igazán szükséges, és a 0.szintet referencia kategóriának kell tekinteni (amely megkönnyíti az eredmények értelmezését):
> könyvtár(reshape2)
> alkohol.df < – melt(alkohol)
> nevek(alkohol.df) < – c (“betegség”,” expozíció”,”n”)
> szintek(alkohol.df$betegség) < – c(1,0)
> alkohol.DF$disease < – relevel(alkohol.df$betegség, “0”)
a logisztikai regressziós modellt ezután a következő módon írják:
> glm.res < – glm(betegség expozíció, adatok = alkohol.DF,
család = binomiális, súlyok = n)
> összefoglaló(glm.res)
az érdeklődés eredménye az expozícióhoz társított sor > = 80 mivel információt nyújt az R által becsült expozíciós változóhoz társított regressziós együttható értékéről, standard hibájával, valamint a vizsgálati statisztika értékéről. Itt a regressziós együtthatót az esélyarány naplójaként értelmezzük. Vegye figyelembe, hogy pontosan ugyanazokat az eredményeket kapnánk, ha kicserélnénk a változók szerepét az előző készítményben, az expozíciós betegségben.
a fent kiszámított esélyarány a kamattényezőhöz (expozíció) kapcsolódó regressziós együtthatóból, valamint annak 95% – os konfidencia intervallumából a következő módon található:
> exp (coef (glm.res))
> exp(confint(glm.res))
a második megoldás az esetek számának és az egyének teljes számának mérlegeléséből áll:
> alkohol2 < – adatok.keret(kiállítások = c(“<> = 80”), case = c(104,96),
a teljes = c(104 + 666, 96 + 109))
> összefoglaló(glm(cbind(az összes eset) kiállítások, adatok = alcohol2,
a családi = binomiális))
4) Ez az adat származik tanulmány célja, hogy létrehozzák a prognosztikai érvényességét kreatin-kináz koncentráció a szervezetben, a megelőzés, az az esemény, a miokardiális infarktus .
Az adatok az sck fájlban érhetők el.dat: az első oszlop megfelel a kreatin-kináz (CK) változónak, a második a változó betegség (pres) jelenlétének, az utolsó pedig a betegség változó hiányának (abs).
A)
mennyi az alanyok teljes száma?
b)
kiszámítja a relatív beteg/egészséges frekvenciákat, és a kreatin-kináz értékek szerinti evolúciójukat scatterplot (pontokat + pontokat összekötő szegmensek) segítségével reprezentálja.
c)
egy logisztikai regressziós modell alapján, amelynek célja a betegség valószínűségének előrejelzése, számítsa ki a CK értékét, amelyből ez a modell azt jósolja, hogy az emberek a betegségben szenvednek, figyelembe véve a 0, 5 küszöbértéket (ha p (beteg) ≥ 0, 5, akkor beteg = 1).
d)
grafikusan reprezentálja a modell által megjósolt rossz valószínűségeket, valamint a CK értékek szerinti empirikus arányokat.
e)
állítsa be a megfelelő ROC görbét, és adja meg a görbe alatti terület értékét.
mivel a változók neve nincs az adatfájlban, az adatok importálása után azonnal hozzá kell rendelni őket:
> sck < – read.táblázat (“sck.dat”, header = FALSE)
> nevek(beteg) < – c(“ck”, “pres”, “abs”)
> összefoglaló(beteg)
A teljes órakeretet megfelel az összeg a frekvenciákat a két változó pres, abs, hogy:
> sum(beteg)
vagy, egyformán: összeg ( sck$pres) + összeg (sck$abs) (de nem sck$pres + sck$abs).
e két változó relatív gyakoriságának kiszámításához meg kell ismerni a teljes mennyiségeket változó szerint. Ezek az apply parancs segítségével érhetők el, és oszlopokkal működtethetők:
> ni < – apply(sck, 2, sum)
ennek alapján elegendő a PRES és abs változók minden értékét a korábban kiszámított értékekkel megosztani. A kapott értékeket két új változóban tároljuk ugyanabban az Adattáblában:
> sck$pres.prop < – sck$pres/ni
> sck$abs.kellék < sck$abs/ni
lehetséges, hogy ellenőrizze, hogy a számításaim helyesek: a értékeinek összege minden változóhoz kell egyenlő 1:
> alkalmazni(beteg, 2, összeg)
Akkor az arányok alapján csupán azt jelenti, hogy egy diagram, figyelembe véve, hogy az értékek a változó ck, mint a abscissae:
> xyplot(pres.prop + abs.prop ck, data = sck, type = c (“b”,”g”),
auto.key = TRUE, ylab = “Frequency”)
az utasítás típusa = c (“b”,” g”) azt jelenti, hogy azt akarjuk, hogy a pontok vonalak (“b” = ” o ” + “l”), valamint egy rács (“g”) szerint jelenjenek meg.
a csoportosított adatok regressziós modellje a következő:
> glm.res < – glm (cbind(pres,abs) ck, data = sck, family = binomial)
> összefoglaló (glm.res)
a valószínűségek formájában kifejezett előrejelzéseket, nem pedig a log-odds skálán, a predict parancs segítségével kapjuk meg a = “válasz” opció típusának megadásával az alábbiak szerint:
> glm.pred < – predict(glm.res, type = “response”)
> names(glm.pred) < – sck$CK
figyelembe véve, hogy a valószínűségek ≥ 0,5 kijelölik a beteg személyeket, így a következő eloszlást kapjuk:
> glm.pred
arra a következtetésre jutottak, hogy az embereket e modell szerint betegnek tekintik a 80 vagy annál nagyobb ck értékekre. Ez könnyen ellenőrizhető egy diagrammal, ahol a CK változó értékei alapján előrejelzett valószínűségek képviseltetik magukat:
> sck$sick < – sck$pres/(sck$pres + sck$abs)
> xyplot(GLM.pred ~ sck$CK, type = “l”,
ylab = “valószínűség”, xlab = “ck”,
panel = function(…) {
panel.xyplot (…)
panel.xyplot (sck$CK, sck$sick, pch = 19, col = “grey”)
})
az első lépésben “dekompresszálni” kell a csoportosított adatokat, és két oszlopból álló táblázatot kell létrehozni: az első a CK változót, a második pedig a betegség jelenlétét vagy hiányát. A számításokat az ni:
> sck változóban korábban kiszámított és elérhető alcsoportonként használjuk.bontsa ki a < – adatokat.keret(ck = c(rep(beteg$ck, beteg$pres),
rep(beteg$ck, beteg$abs)),
a beteg = c(rep(1,ni), képviselő(0,ni)))
> táblázat(beteg.kibontása $ sick)
> a(sck.bontsa ki, tapply(beteg, CK, sum))
az utolsó két parancs célja annak biztosítása, hogy ugyanazok a számok találkozzanak, és hogy a betegek eloszlása a ck értéke szerint helyes legyen. Azt is ellenőrizni, hogy az azonos eredményeket illetően a logisztikai regresszió, ugyanakkor a becsült értékek a korábbi adatok táblázat lehet adni (ugyanaz az eljárás, mint korábban lehet követni a jóslatok lehet reprodukálni, de ez egyszerűbb):
> glm.res2 < – glm(sick CK, data = sck.kibontás, család = binomiális)
> sck.bontsa ki a$prediction < – ifelse(predict(glm.res2,
type = “response”) > = 0.5,
1, 0)
> with(sck.kibontás, táblázat (beteg, jóslat))
az utolsó parancs lehetővé teszi egy olyan zavarmátrix megjelenítését, amelyben a valós diagnosztika és a regressziós modell által előre jelzett értékek keresztbe kerülnek. A helyes osztályozási arányok összehasonlíthatók, ha a referencia küszöbértéket változtatjuk:
> classif.tab < – with(sck.kibontás, táblázat (beteg, jóslat))
> sum(diag(classif.tab)) / sum (classif.tab)
a ROC görbe megjelenítéséhez a rocr csomagot használjuk:
> könyvtár(ROCR)
> pred < – előrejelzés(predict (glm.res2, type = “response”),
sck.expand$sick)
> perf < – performance(pred, “tpr”,”fpr”)
> plot(perf, ylab = “Sensitivity”, xlab = “1-Specificity”)
> grid()
> abline(0, 1, lty = 2)
The value of the area under the curve is obtained as follows:
> performance(pred, “auc”)@”y.values”