8 Tabelanalyser

Vi vil her gennemgå to eksempler med kategoriske data, hvor vi opstiller observationerne i tabeller og finder fordelingen af observationerne med prop.table.

8.1 Envejstabeller

Eksempel: Kønsfordeling blandt børnefødsler (uge 3)

gender = c(rep("D",7), "P")
gender

[1] "D" "D" "D" "D" "D" "D" "D" "P"

Med funktionen table() kan vi få en tabel med antal af hver af de mulige observationer.

tab = table(gender)
tab

gender
D P 
7 1 

Med funktionen prop.table() kan vi beregne andelene af hver type observation:

prop.table(tab)

gender
    D     P 
0.875 0.125 

Drengefødsler udgør således 87.5% af observationerne, og 12.5% af oberservationerne udgøres af pigefødsler.

8.2 CI for andele / binomialtest

Vi undersøger nul-hypotesen:

\(H_0\): Sandsynligheden for drenge- og pigefødsel er lige stor svarende til at der er 50% chance for hhv pige- og dreng.

8.2.1 Eksakt CI + test

Med funktionen binom.test() laves et eksakt test af hypotesen baseret på binomialfordelingen. Her angiver vi antal af hver observation med table(gender) og p = 0.5 angiver at hypotesen, som skal testes, er at sandsynligheden for drenge- hhv. pigefødsel er lige stor. Vi får da givet både p-værdi, estimat af sandsynligheden og konfidensinterval for sandsynligheden.

binom.test( table(gender), p = 0.5)

Hvis vi blot er interesseret i konfidensintervallet, behøver vi ikke specificere p = 0.5 men R udfører alligevel testet per default (dvs helt samme output med binom.test( table(gender) )).

Man kan også bare selv lave optællingen af piger og drenge og angive

binom.test( 7, 8 )

hvor første argument er antal ‘positive/successer’, andet er totalen.

8.2.2 Asymptotisk CI + test

Alternativt kan man bruge funktionen prop.test(), der bygger på en normalfordelings approximation:

prop.test( table(gender), p = 0.5)

prop.test laver pr. automatik en kontinuitetskorrektion. For at undertrykke kontinuitetskorrektionen kan man angive correct = F.

prop.test( table(gender), p = 0.5, correct = F)

Hvis vi blot er interesseret i konfidensintervallet, behøver vi ikke specificere p = 0.5 men R udfører alligevel testet per default (dvs helt samme output med prop.test( table(gender) )).

Man kan også bare selv lave optællingen af piger og drenge og angive

prop.test( 7, 8 )

hvor første argument er antal ‘positive/successer’, andet er totalen.

8.3 Tovejstabeller

Eksempel: Farveblindhed fordelt på køn

Vi observerer forekomsten af farveblindhed blandt 150 drenge og 120 piger. Observationerne opstiller vi i en 2*2-tabel med matrix-funktionen og tilføjer række- og søjlenavne.

fb = matrix( c(1,6,119,144), nrow = 2)
rownames(fb) = c("Piger","Drenge")
colnames(fb) = c("ja","nej")

fb

       ja nej
Piger   1 119
Drenge  6 144

Med prop.table kan vi udregne andelene af farveblindhed bland piger hhv. drenge. Da vores tabel er 2-dimensionel skal vi angive over hvilken margin, der skal beregnes andele. Da vi ønsker at opdele efter køn (rækkerne i tabellen), skal vi angive margin=1. For margin=2 opdeles efter søjler.

prop.table(fb, margin = 1)

             ja    nej
Piger  0.008333 0.9917
Drenge 0.040000 0.9600

Vi kan lave figur, der viser andelene for hvert køn med barplot. For at få den rigtige figur er vi nødt til at vende tabellen først med t-funktionen og så regne andele over søjlerne i stedet.

barplot( prop.table(t(fb), margin = 2) )

For at give overblik over observationerne kan vi udvide tabellen med rækkesummer, søjlesummer samt andelen af piger hhv. drenge med farveblindhed. Dette kan f.eks. gøres vha addmargins().

addmargins(fb)

       ja nej Sum
Piger   1 119 120
Drenge  6 144 150
Sum     7 263 270

8.4 Test for uafhængighed

Eksempel: Farveblindhed

I dette eksempel ønsker vi at undersøge, om risikoen for farveblindhed er uafhængig af køn. Lader vi \(p_1\) og \(p_2\) være risikoen for farveblindhed for hhv. drenge og piger, er vores nul-hypotese:

\(H_0\): \(p_1=p_2\)

Først opstiller vi observationerne i en 2*2-tabel og tilføjer række- og søjlenavne.

fb = matrix( c(1,6,119,144), nrow = 2)
rownames(fb) = c("Piger","Drenge")
colnames(fb) = c("ja","nej")

fb

       ja nej
Piger   1 119
Drenge  6 144

\(\chi^2\)-test for uafhængighed laves med chisq.test. Som i prop.test laver funktionen per automatik en kontinuitetskorrektion. For at undertrykke kontinuitetskorrektionen angives correct = F.

chisq.test(fb)
chisq.test(fb, correct =F)

Hvis tabellen er “tynd” er det bedre at benytte Fishers eksakte test vha. funktionen fisher.test. For at afgøre om vi skal bruge Fishers eksakte test, skal vi beregne de forventede antal observationer i tabellen under antagelse af uafhændighed.

De kan nemt udtrækkes af et chisq.test-objekt ved at bruge $expected. Sidebemærkning: Hvis man vil se en samlet liste over de størrelser, man kan trække ud af et objekt vha. $, kan man bruge funktionen names().

names(chisq.test(fb))

[1] "statistic" "parameter" "p.value"   "method"    "data.name" "observed" 
[7] "expected"  "residuals" "stdres"   

chisq.test(fb)$expected

          ja   nej
Piger  3.111 116.9
Drenge 3.889 146.1

Da der er forventede værdier i tabellen under 5, benytter vi i stedet Eksakt Fisher test med fisher.test, hvorved vi kan aflæse p-værdien.

fisher.test(fb)

    Fisher's Exact Test for Count Data

data:  fb
p-value = 0.1
alternative hypothesis: true odds ratio is not equal to 1
95 percent confidence interval:
 0.004353 1.705674
sample estimates:
odds ratio 
    0.2026 

8.5 Effektmål

Effektmål dækker over:

• Risikodifferens
• Odds-ratio (OR)
• Relativ risiko (RR)

Nedenfor gennemgås et eksempel, hvor effektmålene udregnes og konfidensintervaller bestemmes.

NB: Det er vigtigt at “grupperne” (her piger/drenge) er i rækkerne og “observationerne” (her farveblindhed ja/nej) er i søjlerne. Vær desuden opmærksom på rækkefølgen på søjlerne i tabellen, da det har stor betydning for fortolkningen af effektmålene, der bestemmes nedenfor!

8.5.1 Risikodifferens

Estimat af risikoen for farveblindhed for drenge hhv. piger kan bestemmes med prop.test. Risiko differensen estimeres som differensen mellem de to estimater af sandsynligheder. prop.test giver desuden et konfidensinterval for risiko differensen. Her bestemmes CI for risikodifferensen for at være farveblind (prop.test() baserer beregningerne på første søjle i tabellen:

prop.test(fb)

    2-sample test for equality of proportions with continuity correction

data:  fb
X-squared = 1.5, df = 1, p-value = 0.2
alternative hypothesis: two.sided
95 percent confidence interval:
 -0.07449  0.01116
sample estimates:
  prop 1   prop 2 
0.008333 0.040000 

8.5.2 Relativ risiko og oddsratio

Til at bestemme oddsratio og relativ risiko samt tilhørende konfidensintervaller kan man bruge funktionen epitab() fra pakken epitools. epitools skal først installeres og indlæses:

install.packages("epitools")
library("epitools")

Argumenter til epitab():

Første argument til epitab() er en 2x2 tabel, her fb:

epitab(fb)

Vi kan styre analysen vha ekstra argumenter:

• method angiver metoden, typisk method = "riskratio" eller method = "oddsratio".
• rev som kan benyttes til at bytte rundt på rækker og/eller søjler (REVerse). F.eks. rev=rows, rev=colums eller rev=both
• Vi kan plukke resultaterne ud en tabel ved at benytte $tab, i.e. epitab(fb)$tab

Relativ risiko:

Vi angiver method = "riskratio":

epitab(fb, method = "riskratio")$tab

       ja       p0 nej     p1 riskratio  lower upper p.value
Piger   1 0.008333 119 0.9917    1.0000     NA    NA      NA
Drenge  6 0.040000 144 0.9600    0.9681 0.9333 1.004  0.1364

I søjlerne p0 og p1 får vi rækkeprocenterne. Desuden relativ risiko for ikke-farveblindhed (svarende til nej, med andele i søjle p1) for drenge vs piger til 0.97 (=0.96/0.99). Tilhørende konfidensinterval aflæses under lower og upper.

Vi ønsker i stedet relativ risiko for farveblindhed og kan bytte søjlerne i vores tabel fb med ekstra argument rev="columns":

epitab(fb, method = "riskratio", rev="columns")

$tab
       nej     p0 ja       p1 riskratio  lower upper p.value
Piger  119 0.9917  1 0.008333       1.0     NA    NA      NA
Drenge 144 0.9600  6 0.040000       4.8 0.5858 39.33  0.1364

$measure
[1] "wald"

$conf.level
[1] 0.95

$pvalue
[1] "fisher.exact"

Her vinder vi relativ risiko for farveblindhed for drenge vs piger til 4.8 (=0.04/0.0083)

Oddsratio

Vi ønsker OR for farveblindhed, drenge vs piger, og benytter derfor rev=columns. I søjlerne p0 og p1 får vi rækkeprocenterne.

epitab(fb, method = "oddsratio", rev="columns")$tab

       nej     p0 ja     p1 oddsratio  lower upper p.value
Piger  119 0.4525  1 0.1429     1.000     NA    NA      NA
Drenge 144 0.5475  6 0.8571     4.958 0.5887 41.76  0.1364

Her aflæses OR for farveblindhed for drenge i forhold til piger til 4.96. Tilhørende konfidensinterval aflæses under lower og upper.

8.6 Power

Hvis man kender sandsynligheden for en given hændelse i to forskellige grupper, kan man finde den nødvendige stikprøvestørrelse i hver gruppe for at opnå en ønsket power ved hjælp af power.prop.test.

Eksempel: Farveblindhed blandt drenge og piger

Hvis det antages at risikoen for farveblindhed er 0.01 for piger og 0.04 for drenge, kan man finde den nødvendige stikprøvestørrelse i hver gruppe som nedenfor:

power.prop.test(p1 = 0.01, p2 = 0.04, sig.level = 0.05, power = .80)

     Two-sample comparison of proportions power calculation 

              n = 424
             p1 = 0.01
             p2 = 0.04
      sig.level = 0.05
          power = 0.8
    alternative = two.sided

NOTE: n is number in *each* group

8.7 Parrede binære data

Når vi har parrede binære data, kan vi udføre McNemar test med enten mcnemar.test eller mcnemar.exact.

• mcnemar.test er asymptotisk test
• mcnemar.exact er eksakt test

Eksempel: Sensitivitet af subtrater til detektion af tuberkelbakterier (uge 3)

Resultaterne af dyrkning af de 50 prøver i de to subtrater opstilles i en 2*2 tabel ved først at konstruere en vektor med resultaterne for subtrat a og en vektor med de tilhørende resulter for subtrat b og så anvende funktionen table.

aa = c("+", "+", "-", "-")
bb = c("+", "-", "+", "-")
antal = c(20, 12, 2, 16)

a = rep(aa, antal)
b = rep(bb, antal)
tub = table(a, b)

tub

   b
a    -  +
  - 16  2
  + 12 20

Vi ønsker at undersøge om sensitiviteten (andelen af positive) af de to subtrater er ens vha. et McNemar test.

En asymptotisk McNemar test mcnemar.test kan laves med eller uden kontinuitetskorrektion:

mcnemar.test(tub)
mcnemar.test(tub, correct = F)

For at lave eksakt McNemar test skal vi først installere og indlæse pakken exact2x2.

install.packages("exact2x2")
library(exact2x2)

mcnemar.exact(tub)