2
$\begingroup$

Hi all I have a 2x2 contingency table like this:

            Product B   Not in B    Row Total
Product A       25000       175000      200000
Not in A        4           24899996    24900000
Column Total    25004       25074996    25100000

since I have a cell with less than 5, I have to use the fisher test formula but I have to solve it manually thus I can't used the fisher.test function in R. I have to used this formula http://en.wikipedia.org/wiki/Fisher's_exact_test in order to do the test but factorial function in R and excel can't handle large numbers. Is their a way to compute the fisher exact test without using the factorial or simplifying the equation? Thanks.

Also I'm curious about the equation in fisher.test function in R for it can get the p-value above easily and I can't do it manually due to large factorial. I think that function uses some sort of shortcut. I am still studying it.

$\endgroup$
4
  • 1
    $\begingroup$ It does not matter that one cell has a small count. What matters is its expected value, which is huge. You can trust the chi-squared distribution. $\endgroup$
    – whuber
    Feb 10, 2015 at 3:42
  • $\begingroup$ Thanks @whuber. Forgot that the assumptions only state 80% of the cells must be 5 or more. But I am still looking on how R computes fisher.test. I got a close answer using the algorithms on the code. B=2000 almost.1 <- 1 + 64 * .Machine$double.eps STATISTIC <- sum(lfactorial(x)) (STATISTIC/almost.1)/(B+1) $\endgroup$
    – jbest
    Feb 10, 2015 at 3:54
  • 1
    $\begingroup$ I meant something different: the 80% is irrelevant here, because 100% of your cells have huge expected values. In fact, the departure from independence is so blatantly obvious that any formal test would be superfluous. $\endgroup$
    – whuber
    Feb 10, 2015 at 13:51
  • $\begingroup$ @whuber thank you again. The chi square test and fisher are exactly the same with this kind of data. $\endgroup$
    – jbest
    Feb 11, 2015 at 2:28

1 Answer 1

3
$\begingroup$

The usual rules of thumb about using the chi-square don't refer to the observed data counts, but to the expected counts.

None of your expected counts are small. There's simply no issue here.

As whuber says, you can trust a chi square here -- or you could do a G test -- both will behave similarly in this case. (The Fisher test should also give a similar outcome.)

If you must do a Fisher test, I'd be inclined to use a normal approximation to the hypergeometric; I believe it will give essentially the same result you seek. Alternatively you could simulate to get the p-value.

Note that many packages - R included - have functions to compute log-gamma, log-factorial and log-binomial-coefficient functions. (In R they're lgamma,lfactorial, and lchoose.)


Some comparisons in R:

# set up data
ctb <- matrix(c(175000, 200000, 24899996, 24900000),nr=2,byrow=TRUE)

#chi square
> chisq.test(ctb)

        Pearson's Chi-squared test with Yates' continuity correction

data:  ctb
X-squared = 1654.073, df = 1, p-value < 2.2e-16


 # G test
 expec <- outer(rowSums(ctb),colSums(ctb))/sum(ctb)
 (G <- 2*sum(ctb*log(ctb/expec)))
[1] 1655.443
 pchisq(G,1,lower.tail=FALSE)
[1] 0

# Fisher test. If we proceed naively --
r <- rowSums(ctb)
dhyper(ctb[1,1],r[1],r[2],colSums(ctb)[1])
phyper(ctb[1,1],r[1],r[2],colSums(ctb)[1])

# then these numbers all zeros. This is enough that we can just say p is too small  
# that trying to give an exact value makes no sense. But let's try to compute it 
# anyway by working on the log scale:

dhyper(ctb[1,1],r[1],r[2],colSums(ctb)[1],log=TRUE)
[1] -834.3587

Well, that's a ridiculously low number, but we need it to find what tables are in the "at least as extreme" category. The corresponding smallest value at least as extreme (in terms of likelihood) in the other direction is 199814:

dhyper(199814,r[1],r[2],colSums(ctb)[1],log=TRUE)
[1] -834.4395

So the logs of the required tail probs are:

phyper(ctb[1,1],r[1],r[2],colSums(ctb)[1],log=TRUE)
phyper(199814-1,r[1],r[2],colSums(ctb)[1],log=TRUE,lower.tail=FALSE)

which give -832.2798 and -832.3606 respectively

We'd need to add those two probabilities given on the log-scale. One way to do that is take out a common factor (e.g. shift up by 832 and exponentiate), add the scaled probabilities, take logs and shift back down.

The result gives a combined log-probability of -831.6263, or a p-value of about $10^{-361.17}$. Which we already knew was largely meaningless (for a start, the assumptions won't hold closely enough for that to be even remotely close to a meaningful number).

By contrast the chi-square test gave a log p-value of -830.97 and the unadjusted G test gave a log p-value of -831.65, if I made no errors.

Normal approx to hypergeometric

Mean is nK/N
Variance is nK/N x (N-K)/N x (N-n)/(N-1)

N=sum(ctb)
K=colSums(ctb)[1]
n=rowSums(ctb)[1]
m=n*K/N
v=n*K/N*(N-K)/N* (N-n)/(N-1)
s=sqrt(v)
(z=(ctb[1,1]-m)/s)
pnorm(z,log=TRUE)+log(2) # 2nd term for other tail
-831.0353

With continuity correction:

(z=(ctb[1,1]-m+.5)/s)
pnorm(z,log=TRUE)+log(2) # 2nd term for other tail
-830.97

All these correspond to p-values on the order of $10^{-361}$. It surprises me quite how close together they all are.

[For all of them we might as well just have said "p-value < 2.2e-16" just as the inbuilt chi-square test does.]

$\endgroup$
1
  • $\begingroup$ +1 I'm not in the least surprised all the results are consistently tiny (the p-values actually vary by a factor of four to five), but I am quite surprised they are as large as this! $\endgroup$
    – whuber
    Feb 11, 2015 at 15:06

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.