Vidakovic B. Statistics for Bioengineering Sciences: With Matlab and WinBugs Support

Подождите немного. Документ загружается.

540 14 Models for Tables

The traditional effect size in this context is Cohen’s w, which is in fact

equivalent to the

φ-coefﬁcient, φ =

/n, where n is the total table size. Un-

like the

φ-coefﬁcient which is used for 2×2 tables only, w is used for arbitrary

×c tables and can exceed 1. Effects w =0.1,0.3 and 0.5 correspond to a small,

medium and large size, respectively.

For prospective analyses, the noncentrality parameter

λ is nw

while for

retrospective analyses

λ is the observed χ

. The power is

−β =1 −ncχ

(χ

k,1

−α

, k,λ),

where k

=(r −1)×(c −1) is the number of degrees of freedom and χ

k,1

−α

is the

−α)-quantile of a (central) χ

distribution.

Example 14.5. (a) Find the power in a 2

×6 contingency table, for n = 180 and

=0.3 (medium effect).

w = 0.3; n = 180; k = (2-1)

(6-1); alpha = 0.05; lambda=n

w^2

pow = 1-ncx2cdf( chi2inv(1-alpha,k), k, lambda) %0.8945

(b) What sample size ensures the power of 95% in a contingency table 2×6,

for an effect of w

=0.3 and α =0.05?

beta = 0.05; alpha = 0.05; k = (2-1)

(6-1); w=0.3;

pf = @(n) ncx2cdf( chi2inv(1-alpha,k), k, n

w^2) - beta;

ssize = fzero(pf, 200) %219.7793 approx 220

14.2.2 Cohen’s Kappa

Cohen’s kappa is an important descriptor of agreement between two testing

procedures. This descriptor is motivated by calibrating the observed agree-

ment by an agreement due to chance. If p

is the proportion of agreement due

to chance and p

the proportion of observed agreement, then

κ =

− p

1 − p

. (14.2)

For a paired table representing the results of n tests by two devices or

ratings by two raters

+ −

+ a b a +b

− c d c +d

a +c b +d n = a +b +c +d

Cohen’s kappa index is deﬁned as

14.2 Contingency Tables: Testing for Independence 541

κ =

2(ad −bc)

(a +b)(b +d) +(a +c)(c +d)

. (14.3)

The expression in Eq. (14.2) is equivalent to that in Eq. (14.3) for p

=a/n+

d/n and p

=(a +b)(a +c)/n

+(b +d)(c +d)/n

, respectively. The former is the

observed agreement equal to the proportion of (+,+) and (–,–) outcomes, while

the latter is the proportion of agreement when the results are independent

within ﬁxed marginal proportions.

There are no formal rules for judging

κ, but here is the standard:

κ Degree of agreement

<0.20 Poor

0.20 – 0.40 Fair

0.40 – 0.60 Moderate

0.60 – 0.80 Good

0.8 – 1 Very good

The MLE of κ is

mle

4(ad −bc) −(b −c)

(2a +b +c)(2d +b +c)

and it is obtained from (14.2) by taking p

= a/n + d/n and p

= P

+(P

)

where P

=(2a+b+c)/(2n) and P

=(2d +b+c)/(2n) are, respectively, the MLEs

of the prevalence of + and – in the population.

There are several expressions for the variance of

κ. Two standard estima-

tors are the Block–Kraemer (BK) and Garner (G) approximations:

(1) (BK):

Var

κ ≈

1 −

(1 −

κ)(1 −2

κ) +

κ(2 −

κ)

2P(1 −P

)

;

(2) (G)

Var

κ ≈

(1 − p

)

(

a+1

b+1

c+1

d+1

)

The sampling distribution of

κ is asymptotically normal and the approx-

imation is satisfactory if n is large and

κ not too close to 1. This leads, in a

standard manner, to approximate conﬁdence intervals for

κ.

Example 14.6. A company producing a medical sensor A is applying for FDA

approval of a new version B. Both sensors A and B are prone to errors, and a

542 14 Models for Tables

gold standard is absent. The FDA is requesting that the new sensor be compa-

rable to the currently used one, and the company decides to include Cohen’s

statistic in the report.

The experiment consisted of n

=2803 trials and resulted in a paired table

[97 11; 6 2689], where 97 was the number of (+,+) outcomes and 2689 the

number of (–,–) outcomes. The code

cohen.m ﬁnds Cohen’s

κ and 95% conﬁ-

dence intervals for the population

κ, based on the two estimators of variance,

Block-Kraemer and Garner.

data =[97 11; 6 2689]

%data =

% 97 11

% 6 2689

a=data(1,1); b=data(1,2); c=data(2,1); d=data(2,2);

apb = a+b; cpd = c+d; apc = a+c; bpd = b+d;

n = a + b + c + d;

%------------

p0 = (a + d)/n; %Observed agreement

pc = (apb

apc + cpd

bpd)/n^2; %Chance agreement

Pre = (2

a + b + c)/(2

n); %Prevalence of +

pcc = Pre^2 + (1-Pre)^2; %MLE of chance agreement

%------------

kappa = 2

d - b

c)/(apb

bpd + apc

cpd) %0.9163

%or kappa=(p0-pc)/(1-pc)

kappamle = (4

d - b

c) - (b -c)^2)/...

((2

a + b + c)

d + b + c)) %0.9163

%or kappamle=(p0-pcc)/(1-pcc)

%------------

%Block-Kraemer variance estimator

varbk = (1- kappa)/n

( (1- kappa)

(1-2

kappa) + ...

(kappa

(2 - kappa)/(2

Pre

(1-Pre)))) %4.0731e-004

%Garner variance estimator

vargarner = 4/( (1- pc)^2

n^2

(1/(a+1) + 1/(b+1) + ...

1/(c+1) + 1/(d + 1) ) ) %4.0971e-004

%------------

%Confidence intervals

[ kappa - 1.96

sqrt(varbk) ...

kappa + 1.96

sqrt(varbk)] %0.8767 0.9558

[ kappa - 1.96

sqrt(vargarner) ...

kappa + 1.96

sqrt(vargarner)] %0.8766 0.9560

Cohen’s κ is estimated to be 91.63%, which represents very good agree-

ment.

14.3 Three-Way Tables 543

14.3 Three-Way Tables

A natural extension of two-way tables for testing the independence of two fac-

tors are n-dimensional tables for testing the independence of n factors. We

will discuss a three-dimensional extension; the interested reader can consult

Zar (2007), Agresti (2002), or Fienberg (2000) for more detailed coverage. In

three-dimensional tables the counts constitute three-dimensional arrays char-

acterized by rows, columns, and pages. We associate factors with these three

dimensions and consider row, column, and page factors (R, C, and P). In the

cell (i, j, k) there are n

i jk

observations, i =1,... , r; j =1, . .., c; and k =1, . .., p.

Denote the total number of observations by n. The empirical probabil-

ity of the cell (i, j, k) is n

i jk

/n, and the empirical marginal probabilities of

row i, column j, and page k are n

i··

/n, n

·j·

/n, and n

··k

/n. The numerators

are calculated as the sums over all indices replaced by dots. For example,

·j·

i jk

, j =1,.. . , c.

We are interested in testing the hypothesis H

that the factors R, C, and P

are independent. The alternative H

would be that the factors are dependent.

Under H

, the frequency in the cell (i, j, k) is expected to be

i jk

= n ×

i··

·j·

··k

as the product of the total number of observations n and the corresponding

empirical marginal probabilities. After simpliﬁcation, the expected frequency

in the cell (i, j, k) becomes

i jk

i··

×n

·j·

×n

··k

The test statistic is

i=1

j=1

k=1

i jk

−e

i jk

)

i jk

, (14.4)

and the likelihood ratio statistic is

i=1

j=1

k=1

i jk

log

i jk

544 14 Models for Tables

When H

is true, both statistics χ

and G

follow approximately a χ

distribu-

tion with rc p

−r − c − p +2 degrees of freedom. Large values of χ

are critical

for H

If H

is rejected, then multiple alternatives are possible. For example, all

three factors R, C, and P are mutually dependent, factors R and C are de-

pendent but both are independent of factor P, and so on. To specify why H

rejected, one needs to test the hypothesis whether each single factor is inde-

pendent of the other two. There are three such tests, and they are summarized

in the following table:

Factor e

i jk

d f

R vs. C, P

×n

·jk

rc p −cp −r +1

C vs. R, P

×n

·jk

rc p −r p −c +1

P vs. R, C

×n

·jk

rc p −rc − p +1

Here n

i j·

i jk

, n

·jk

i jk

, and n

i·k

i jk

. Also, as before n

i··

i jk

, n

·j·

i jk

, and n

··k

i jk

. The χ

statistic is calculated

as in (14.4) and the degrees of freedom are given in the table.

Example 14.7. Anolis Lizards of Bimini. This well-known data set comes

from the paper of Schoener (1968) and is also used in Fienberg (1970). The

researcher was interested in structural habitat categories for Anolis lizards

of Bimini: sagrei (brown anole) adult males versus distichus (trunk anole)

adult and subadult males. The brown anole and trunk anole are medium-

sized, fairly robust, “trunk-ground” lizards. They generally prefer the fairly

open vegetation of disturbed sites, where they adopt a head-down, sit-and-

wait posture and perch low on large trunks or fenceposts (Fig. 14.2).

Fig. 14.2 Brown anole (Anolis sagrei) in typical position at its perch.

The researcher was interested in the preferences of these two species with

respect to the perch height and diameter.

14.3 Three-Way Tables 545

A. sagrei A. distichus

Perch Diameter Perch Diameter

≤4 >4 ≤4 >4

Perch Height >4.75 32 11 61 41

(in feet)

≤4.75 86 35 73 70

The data is classiﬁed with respect to three dichotomous factors: Height at

levels Low (

≤4.75) and High (>4.75), Diameter at levels Small (≤4) and Large

(

> 4), and Species with levels A. sagrei and A. distichus. Are these three fac-

tors independent? The null hypothesis is that the factors are independent and

the alternative is that they are not. The MATLAB program

tablerxcxp.m

calculates the expected frequencies (under mutual independence hypothesis)

and provides the

statistic, its degrees of freedom, and the p-value. This

program also outputs the expected frequencies in the format that matches the

input data. Here, Height is the row factor R, Diameter is the column factor C,

and Species is the page factor P.

anolis = [32 11; 86 35];

anolis(:,:,2)=[61 41; 73 70];

[ch2 df pv exp]=tablerxcxp(anolis)

%ch2 = 23.9055

%df = 4

%pv =8.3434e-005

%exp(:,:,1) =

% 35.8233 22.3185

% 65.2231 40.6351

%exp(:,:,2) =

% 53.5165 33.3417

% 97.4370 60.7048

The hypothesis H

is rejected (with p-value < 5%), and we infer that the

three factors are not independent. However, this analysis does not fully ex-

plain the dependencies responsible for rejecting H

. Are all three factors de-

pendent, or maybe two of the factors are mutually dependent and the third

is independent of both? When H

is rejected, we need a partial independence

test, similar to pairwise comparisons in the case where the ANOVA hypothesis

is rejected. The MATLAB program for this test is

partialrxcxp.m.

[ch2 df pv exp]=partialrxcxp(anolis,’r

cp’)

%ch2 = 12.3028; df=3; pv = 0.0064, exp=...

[ch2 df pv exp]=partialrxcxp(anolis,’c

pr’)

%ch2 = 14.4498; df=3; pv = 0.0024, exp=...

[ch2 df pv exp]=partialrxcxp(anolis,’p

rc’)

%ch2 = 23.9792; df=3; pv =2.5231e-005, exp=...

Three tests are performed: (i) the row factor independent of column/page

factors (

partial = ’r

cp’), (ii) the column factor independent of row/page fac-

tors (

partial = ’c

rp’), and (iii) the page factor independent of row/column

factors (

partial = ’p

rc’). It is evident from the output that all three tests

546 14 Models for Tables

produced signiﬁcant χ

, that is, each factor depends on the two others. Since

Height was the row factor R, Diameter the column factor C, and Species the

page factor P, we conclude that the strongest dependence is that of Species

factor on the height and diameter of the perch,

partialrxcxp(anolis,’p

rc’),

with a p-value of 0.000025.

14.4 Contingency Tables with Fixed Marginals: Fisher’s

Exact Test

In his book, Fisher (1935) provides an example of a small 2 ×2 contingency

table related to a tea-tasting experiment, namely, a woman claimed to be able

to judge whether tea or milk was poured in a cup ﬁrst. The woman was given

eight cups of tea, in four of which tea was poured ﬁrst, and was told to guess

which four had tea poured ﬁrst. The contingency table for this design is

Guess milk ﬁrst Guess tea ﬁrst Total

Milk ﬁrst x 4−x 4

Tea ﬁrst 4 −x x 4

Column total 4 4 8

The number of correct guesses “Milk ﬁrst” in the cell (1,1), x, can take

values 0, 1, 2, 3, or 4 with the probabilities

(

)(

−x

)

(

)

, as in

hygepdf(0:4, 8,4,4)

% ans = 0.0143 0.2286 0.5143 0.2286 0.0143

These are hypergeometric probabilities, applicable here since the marginal

counts are ﬁxed.

For x

=4 the probability of obtaining this table by chance is 0.0143, while

for

3 the probability of getting this or a more extreme table by chance

is 0.2286 + 0.0143 = 0.2429, and so on. Thus the probability of the woman’s

guessing correctly, i.e., if x

=4, would be less than 5%.

Suppose that in the table

Column 1 Column 2 Row total

Row 1 a b a +b

Row 2 c d c +d

Column total a +c b +d n = a +b +c +d

marginal counts a +b, c +d, a +c, and b +d are ﬁxed. Then a, b, c, and d are

constrained by these marginals. We are interested if the probabilities that an

observation will be in column 1 are the same for rows 1 and 2. Denote these

probabilities as p

and p

. The null hypothesis here is not the hypothesis of

independence but the hypothesis of homogeneity, H

: p

= p

. The test is close

14.4 Fisher’s Exact Test 547

to the two-sample problem considered in Chap. 10; however, in this case the

samples are dependent because of marginal constraints.

The statistic T to test H

is simply the number of observations in the cell

(1,1):

T =a.

If H

is true, then T has a hypergeometric distribution H G (n, a+c,a+b), that

P(T = x) =

a+c

b+d

+b−x

, x =0,1,..., min{a +b, a +c}.

The p-value against the one-sided alternative H

: p

< p

is hygecdf(a,

n, a+c, a+b)

and against the alternative H

: p

> p

is 1 - hygecdf(a-1, n,

a+c, a+b)

. If the hypothesis is two-sided, then the p-value cannot be obtained

by doubling one-sided p-value due to asymmetry of the hypergeometric dis-

tribution. The two-sided p-value is obtained as the sum of all probabilities

hygepdf(x, n, a+c, a+b), x = 0,1, ...,min{a +b, a + c} which are smaller than

or equal to

hygepdf(a, n, a+c, a+b).

Example 14.8. There are 22 subjects enrolled in a clinical trial and 9 are fe-

males. Researchers plan to administer 11 portions of a drug and 11 placebos.

Only 2 females are administered the drug. Are the proportions of males and

females assigned to the drug signiﬁcantly different? What are the p-values for

one- and two-sided alternatives?

a = 2; b = 7; c = 9; d =4;

n = a + b + c + d;

T = a;

pval = hygecdf(T,n,a+c,a+b) %H1: p1<p2 pval=0.0402

pa = hygepdf(T,n,a+c,a+b);

for i = 1:min(a+b, a+c)+1

p(i) = hygepdf(i-1,n,a+c,a+b) ;

end

pval2 = sum(p(p <= pa)) %H1: p1 ~= p2 pval2=0.1179

Since the one-sided p-value is less than 5%, we reject the hypothesis of

homogeneity of adminstration of a drug versus placebo with respect to gender.

For the two-sided alternative, we fail to reject H

Example 14.9. The Effect of Passive Smoking on Lung Cancer. Lawal

(2003) considers the following data originally published by Correa et al. (1983)

548 14 Models for Tables

on the effect of passive smoking on lung cancer. A total of 155 non-smoking

ever-married females were tabulated by their lung cancer status and hus-

band’s smoking status.

Is the proportion of lung cancer cases homogeneous with respect to the

husband’s smoking status? Find the p-value for both one- and two-sided alter-

natives.

Smoking status

Case Control Total

Spouse smoked 14 61 75

Spouse did not smoke 8 72 80

Total 22 133 155

Here H

: p

= p

and H

: p

> p

or H

: p

6= p

a = 14; b = 61; c = 8; d =72;

n = a + b + c + d;

T = a;

%H1: p1 > p2

pval = 1-hygecdf(T-1,n,a+c,a+b) %0.0941

%H1: p1 ~= p2

pa = hygepdf(T,n,a+c,a+b);

for i = 1:min(a+b, a+c)+1

p(i) = hygepdf(i-1,n,a+c,a+b) ;

end

pval2 = sum(p(p <= pa)) %0.1669

Thus, Fisher’s exact test fails to reject the null hypothesis at 5% signiﬁ-

cance level.

Fisher’s exact test remains valid for designs with random row totals, ran-

dom column totals, or tables with random marginals (as in the previous sec-

tion). In this case the tests are conservative, and more powerful versions exist.

A beneﬁt of using Fisher’s exact test is that it operates with small cell frequen-

cies, for example, a 0 count in a table cell is a possibility. Since

or normal

approximations assume large n and np

i j

s preferably larger than 5, the reason

for the popularity of Fisher’s exact test is obvious.

14.5 Multiple Tables: Mantel–Haenszel Test

In Chap. 10, p. 380, 2×2 tables were discussed in the context of comparing two

proportions. Here we discuss multiple 2

×2 tables and inference from combined

information. The Mantel–Haenszel (Fig. 14.3a,b) methodology can be used in

×2 tables to control for a variable that stratiﬁes the data. This leads to mul-

tiple tables, one for each level of controlled variable. The Mantel–Haenszel

methodology can be used for (i) testing the conditional independence of two

14.5 Multiple Tables: Mantel–Haenszel Test 549

factors or (ii) measuring the degree of conditional association (risk ratios).

All conditioning is on the variable by which the tables are stratiﬁed. There

are several other uses of the Mantel–Haenszel methodology such as in sur-

vival analysis (longrank tests of Peto and Peto) and in depairing of McNemar’s

paired designs.

14.5.1 Testing Conditional Independence or Homogeneity

Suppose that k independent classiﬁcations into a 2 ×2 table are observed. We

could denote the ith such table by

The tables give counts broken down by binary levels of two factors, and the

separate tables usually correspond to the levels of a third factor that needs

to be controlled. Imagine that we want to test for the independence of two

factors, say, political association (Democrat, Republican) and opinion about

some social issue (Support, Oppose). The single contingency table may not be

signiﬁcant, but when controlled by gender (two tables, one for males, the other

for females) or by age group, the dependence may turn out signiﬁcant. Thus,

multiple tables make inference more precise by controlling for an inﬂuential

variable.

For each of k tables consider cell at the position (1,1), so called pivot, with

counts. If the two tabulated factors are independent, then the counts a

should be close to “expected” counts e

=(a

)(a

)/n

The test statistic measuring discrepancies in the pivot cell over all k tables

(|A −E|−1/2)

, where (14.5)

i=1

, E =

i=1

, V =

i=1

)(c

)(a

)(b

)

−1)

which has an approximately

-distribution with 1 degree of freedom when

the (null) hypothesis of independence/homogeneity is true. Large values of

are critical for H