turtlegraphics - 11 months ago 40

R Question

I’m confused about the operation of binom.test.

Say I want to test a sample of 4/10 success against p=0.5.

The P value should be:

P(X <= 4) + P(X >=6)

or

P(X <= 4) + 1-P(X <= 5)

and indeed:

`>pbinom(4,10,p=0.5) + 1-pbinom(5,10,0.5)`

[1] 0.7539063

or:

`>binom.test(4,10,p=0.5)`

Exact binomial test

data: 4 and 10

number of successes = 4, number of trials = 10, p-value = 0.7539

But now I want to test a sample of 95/150 against p=0.66

Here, the expected value is 99, so the P value should be

P(X <= 95) + P(X >= 103)

or

P(X <= 95) + 1-P(X <= 102)

which is

`>pbinom(95,150,.66) + 1-pbinom(102,150,.66)`

[1] 0.5464849

but

`>binom.test(95,150,.66)`

Exact binomial test

data: 95 and 150

number of successes = 95, number of trials = 150, p-value = 0.4914

In fact, the difference in the two P-values is exactly

`dbinom(103,150,.66)`

The only explanation I can guess for this is that there is roundoff error due to the inexact representation of .66 causing R to just miss X=103. Is this all it is, or is there something else going on?

Answer Source

Here is the code for computing the p-value in binom.test(x = 95, n = 150, p= 0.66)

```
relErr <- 1 + 1e-07
d <- dbinom(x, n, p)
m <- n * p
i <- seq.int(from = ceiling(m), to = n)
y <- sum(dbinom(i, n, p) <= d * relErr)
pbinom(x, n, p) + pbinom(n - y, n, p, lower.tail = FALSE)
```

So, binom.test doesn't look symmetrically about the expected value. It looks for the first integer C such that C is bigger than or equal to the expected value and the probability of getting exactly C successes is less than or equal to the probability of getting exactly x successes, up to the fudge factor in relErr. So, instead of saying that p is the probability of getting "at least that far away from the expected value", they say that p is the probability that the probability is at least as small as the value that you obtained.

In this case,

```
dbinom(95,n,p)
```

is 0.05334916. So, binom.test looks for the values of x such that dbinom(x,n,p) is less than 0.05334916. It turns out that those are 0:95 and 104:150. So, binom.test returns the value of

```
sum(dbinom(0:95,n,p)) + sum(dbinom(104:150,n,p))
```

which is 0.4914044.