11
$\begingroup$

I'd like to do a power analysis for a single sample from binomial data, with $H_0: p = 0$, vs. $H_1: p = 0.001$, where $p$ is the proportion of successes in the population. If $0 < p <1$, I could use either the normal approximation to binomial, or $\chi^2$-test, but with $p =0$, these both fail. I'd love to know if there is a way to do this analysis. I'd very much appreciate any suggestions, comments, or references. Many thanks!

$\endgroup$
2
  • $\begingroup$ So why do you not use the exact Clopper-Pearson test ? $\endgroup$ Commented Oct 2, 2012 at 5:22
  • 3
    $\begingroup$ I hope you have a really big sample! This is going to be hard to test. $\endgroup$
    – Peter Flom
    Commented Oct 2, 2012 at 10:32

3 Answers 3

13
$\begingroup$

You have a one-sided, exact alternative hypothesis $p_{1} > p_{0}$ where $p_{1} = 0.001$ and $p_{0} = 0$.

  • The first step is to identify a threshold $c$ for the number of successes such that the probability to get at least $c$ successes in a sample of size $n$ is very low under the null hypothesis (conventionally $\alpha = 0.05$). In your case, $c=1$, regardless of your particular choice for $n \geqslant 1$ and $\alpha > 0$.
  • The second step is to find out the probability to get at least $c$ successes in a sample of size $n$ under the alternative hypothesis - this is your power. Here, you need a fixed $n$ such that the Binomial distribution $\mathcal{B}(n, p_{1})$ is fully specified.

The second step in R with $n = 500$:

> n  <- 500                 # sample size
> p1 <- 0.001               # success probability under alternative hypothesis
> cc <- 1                   # threshold
> sum(dbinom(cc:n, n, p1))  # power: probability for cc or more successes given p1
[1] 0.3936211

To get an idea how the power changes with sample size, you can draw a power function: enter image description here

nn   <- 10:2000                 # sample sizes
pow  <- 1-pbinom(cc-1, nn, p1)  # corresponding power
tStr <- expression(paste("Power for ", X>0, " given ", p[1]==0.001))
plot(nn, pow, type="l", xaxs="i", xlab="sample size", ylab="power",
     lwd=2, col="blue", main=tStr, cex.lab=1.4, cex.main=1.4)

If you want to know what sample size you need to achieve at least a pre-specified power, you can use the power values calculated above. Say you want a power of at least $0.5$.

> powMin <- 0.5
> idx    <- which.min(abs(pow-powMin))  # index for value closest to 0.5
> nn[idx]     # sample size for that index
[1] 693

> pow[idx]    # power for that sample size
[1] 0.5000998

So you need a sample size of at least $693$ to achive a power of $0.5$.

$\endgroup$
10
  • $\begingroup$ According to pwr.p.test, for a power of 0.5, you need at least 677 observations. But power = 0.5 is very low! $\endgroup$
    – Jessica
    Commented Oct 2, 2012 at 14:09
  • $\begingroup$ @caracal Are you using a normal approximation to get your power curve? An exact binomial power function would not be so smooth. It is actually sawtoothed which you could see if the sample size axis is magnified. I discuss this in my 2002 paper in the American Statistician coauthored with Christine Liu. Also the binomial is so highly skewed at very low p that n must be large for the normal approximation to work well. $\endgroup$ Commented Oct 2, 2012 at 15:30
  • 2
    $\begingroup$ @MichaelChernick No, this is from the binomial distributions, not from a normal approximation. Of course you're right that - in general - power for a binomial test is a sawtooth function that is not monotonic. But note that we have a special case here with $p_{0} = 0$. This means that the acceptance region for the alternative hypothesis always starts at 1, regardless of $n$. With a constant threshold $c = 1$, a constant $p_{1} = 0.001$, power is a strictly inreasing function of $n$. $\endgroup$
    – caracal
    Commented Oct 2, 2012 at 19:04
  • $\begingroup$ @Jessica Note that pwr.p.test() uses a normal approximation, not the exact binomial distributions. Just type pwr.p.test to have a look at the source code. You'll find the calls to pnorm() indicating that an approximation is used. $\endgroup$
    – caracal
    Commented Oct 2, 2012 at 19:26
  • 1
    $\begingroup$ @caracal So can I look at it this way: Under the null hypothesis the probability of success is 0 hence if you ever see a success you can reject the null hypothesis. That is why you say the threshold is 1 because if the binomial sum ever gets to 1 you can reject with a type 2 error of 0! Now under the alternative the probability of the first success on the nth trial is (1-p)$^n$$^-$$^1$ p. This probability goes to 0 as n goes to infinity. So a sequential rule going by stop when S$_n$=1 would have power 1 for any p>0. $\endgroup$ Commented Oct 2, 2012 at 19:27
3
$\begingroup$

You can answer this question easily with the pwr package in R.

You will need to define a significance level, power, and effect size. Typically, significance level is set to 0.05 and power is set to 0.8. Higher power will require more observations. Lower significance level will decrease power.

The effect size for proportions used in this package is Cohen's h. The cutoff for a small h is often taken to be 0.20. The actual cutoff varies by application, and might be smaller in your case. Smaller h means more observations will be required. You said your alternative is $p = 0.001$. That is very small

> ES.h(.001, 0)
[1] 0.0632561

But we can still proceed.

 > pwr.p.test(sig.level=0.05, power=.8, h = ES.h(.001, 0), alt="greater", n = NULL)

 proportion power calculation for binomial distribution (arcsine transformation) 

          h = 0.0632561
          n = 1545.124
  sig.level = 0.05
      power = 0.8
alternative = greater

Using these values, you need at least 1546 observations.

$\endgroup$
0
1
$\begingroup$

In your specific case there is a simple exact solution:

Under the particular null hypothesis $H_0: p=0$ you should never observe a success. So as soon as you observe one success you can be sure that $p\neq0$.

Under the alternative $H_1: p=0.001$ The number of trials required to observe at least 1 success follows a geometric distribution. So in order to obtain the minimum sample size to achieve a power of $1-\beta$, you need to find the smallest k such that, $$1-\beta \leq 1-(1-p)^{(k-1)}$$

So with $p=0.001$ to get $80%$ power you would need at least 1610 samples.

$\endgroup$
1
  • $\begingroup$ On reading the comments to solution 1, I realize that this is essentially the same solution as the one you get if sticking to answere one. Nevertheless, it never harms to spell out some basic probability theory results, without the need of arriving there by intuition. $\endgroup$
    – float
    Commented Jun 29, 2017 at 9:10

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.