Re: [R] power 2x3 exact test

From: Ted Harding <ted.harding_at_nessie.mcc.ac.uk>
Date: Thu, 10 May 2007 00:12:48 +0100 (BST)


On 09-May-07 22:00:27, Duncan Murdoch wrote:
> On 09/05/2007 5:11 PM, Bingshan Li wrote:
> > Hi, all,
> >
> > I am wondering if there is an algorithm for calculating power of 2x3
> > table using Fisher exact test. There is one (algorithm 280) for 2x2
> > Fisher exact test but I couldn't find one for 2x3 table. If we are
> > not lucky enough to have one, is there any other way to calculate
> > exact power of 2x3 table? The reason why I want exact power is
> > because some cells are assumed to be very small and chi square
> > approximation is not valid.
>
> I think there are lots of possible alternatives to the null in a 2x3
> table, so you may have trouble finding a single answer to this
> question.
> But assuming you have one in mind, I'd suggest doing a Monte Carlo
> power calculation: simulate a few thousand tables from the alternative
> distribution, and see what the distribution of p-values looks like.
>
> Duncan Murdoch

I'd back Duncan on that point!

More specifically, for the 2x2 table, the table, conditional on the marginals, is a function of one element (say top left-hand corner), and the probability of any table depends on the single parameter which is the odds-ratio of the 4 cell probabilities.

so this case is relatively easy and straightforward and, indeed, for the 2x2 table R'a fisher.test() allows you to specify the odds-ratio as a "null" parameter.

This is not the case with fisher.test() for a larger (say 2x3 table), so to investigate that case you cannot use fisher.test().

In all cases, however (according to the FORTRAN code on which itis based -- see the reference in "?fisher.table"), the rejection region for the exact fisher.test() consists of those table with the smallest probabilities.

For the 2x3 table, say (cell counts with margins, and probabilities):

   a1  b1  c1 | d1          p1  q1  r1
              |
   a2  b2  c2 | d2          p2  q2  r2


   a b c | n

so that

   a1+b1+c1 = d1, a2+b2+c2 = d2,
   a1+a2 = a, b1+b2 = b, c1+c2 = c

the table is a function of any two functionally independent cells (say a1 and b1), and its probability is a function of some two odds-ratios, say

   (p1*r2)/(r1*p2)

   (q1*r2)/(r1*q2)

which, for the standard null hypothesis, are both equal to 1. However, as Duncan says, alternatives are 2-dimensional and so there is not a unique natural form for an alternative (as opposed to the 2x2 case, where it boils down to (p1*q2)/(p2*q1) being not equal to 1, therefore greater than 1, or less than 1, or 2-sidedly either >1 or <1).

The probability of the 2x3 table is proportional to

  ((p1*r2)/(r1*p2))^a1 * ((q1*r2)/(r1*q2))^b1

(or equivalent), divided by the product of the factorials of a1, b1, c1, a2, b2, c2, subject to summing to 1 over all combinations of (a1,b1) giving rise to a table compatible with the marginal constraints.

Given that you expect some cells to be small, it should not be a severe task to draw up a list of (a1,b1) values which correspond to rejection of the null hypothesis (that both ORs equal 1), and then the simulation using different values of the two odds-ratios will give you the power for each such pair of odds-ratios.

The main technical difficulty will be simulation of random tables, conditional on the marginals, with the probabilities as given above.

I don't know of a good suggestion for this. The fisher.test() function will not help (see above). In the case of the 2x2 table, is is a straightforward hypergeometric distribution, but 2x3 is not straightforward. Maybe someone has written a function for this kind of application, and can point it out to us???

A quick R-site search did not help!

Best wishes,
ted.



E-Mail: (Ted Harding) <Ted.Harding_at_manchester.ac.uk> Fax-to-email: +44 (0)870 094 0861
Date: 10-May-07                                       Time: 00:12:29
------------------------------ XFMail ------------------------------

______________________________________________
R-help_at_stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Wed 09 May 2007 - 23:19:48 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Thu 10 May 2007 - 16:31:50 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.