1

I am using prop.test() to calculate the confidence interval of a proportion. I found that prop.test() result differs from the z-score confidence interval when $\hat p$ is small, but coincide otherwise. In particular, with $\hat p=0.916$ it works well:

n = 500
p_hat = 0.916
alpha = 0.1
se = sqrt(p_hat*(1-p_hat)/n)
z = -qnorm(alpha/2)
IC = c(p_hat-z*se,p_hat+z*se); IC
## [1] 0.8955953 0.9364047

prop.test(p_hat*n,n=n,alternative="two.sided", conf.level=1-alpha, correct=FALSE)
## 0.8932886 0.9342336

But with $\hat p = 0.016$ the results differ:

IC = c(p_hat-z*se,p_hat+z*se); IC
## [1] 0.006770041 0.025229959

prop.test(p_hat*n,n=n,alternative="two.sided", conf.level=1-alpha, correct=FALSE)  
## 0.009038314 0.028171427

Do you know how in particular prop.test() calculates the confidence interval?

BruceET
  • 52,418
Freeman
  • 201
  • The CI given by prop.test inverts the test to give a generally more accurate CI than the 'Wald score' CI $\hat p \pm 1.96\sqrt{\hat p(1-\hat p)/n},$ in which SE is estimated by using $\hat p$ instead of $p.$ See documentation from R console window at ? prop.test. // The Wald interval is asymptotically correct as promised, but can fail badly even for moderately large $n.$ Difficulties are that it uses a normal approximation to binomial along with an estimated SE. (Google 'Agresti' for more accurate 95% CI.) – BruceET Nov 18 '17 at 05:02

1 Answers1

1

With $X = 484$ successes in $n = 500$ binomial trials, here is a comparison of the Agresti-style 95% CI with the output of 'prop.test`. The so-called Wilson CI may work even better.

n = 500;  x = 484; p0m = -1:1
x.a = x+2; n.a = n+4  # Agresti: append 2 S's & 2 F's
p.a = x.a/n.a;  se.a = sqrt(p.a*(1-p.a)/n.a)
p.a + p0m*1.96*se.a
## 0.9480839 0.9642857 0.9804876

prop.test(484, 500, correct=FALSE)$conf
## 0.9486551 0.9802085
## attr(,"conf.level")
## 0.95

Notes: (1) Getting good CIs by any method can be challenging when the proportion of successes is very near 0 or 1.

(2) The Agresti CI is not just an ad hoc trick. By conflating $1.96$ with $2$s that occur when inverting the test (and ignoring a tiny term) a 95% Agresti CI comes very close to the test-inverted CI, and yet is simple enough to remember.

(3) A frequentist re-purposing of a Bayesian probability interval based on a $\mathsf{Unif}(0,1) \equiv \mathsf{Beta}(1,1)$ prior, gives the 95% CI:

qbeta(c(.025,.975), x+1, n-x+1)
## 0.9486549 0.9801114

(4) Perhaps also see this page.

BruceET
  • 52,418
  • Thank you for your clarification regarding the improved performance of Agresti's method. I didn't know it.

    However, your calculated CIs: ## 0.9480839 0.9642857 0.9804876, ## 0.9486549 0.9801114, both differ from the CI reported by prop.test: ## 0.9486551 0.9802085.

    Do you know what calculations prop-test() does to compute the CI?

     – Freeman Nov 19 '17 at 09:46
  • Differences beyond the 3rd decimal place can be due to rounding errors. The prop.test interval is from 'inverting' the test. You can do a grid search to verify that to several places. Agresti is close to that but not exactly. You can google 'binomial confidence interval' for Agresti, Wilson, etc. methods. Endless discussion about this in stat literature. Once again, all methods have serious weaknesses (computational and philosophical) when the pop prop of successes in very near 0 or 1. Some intervals can stray outside $(0,1)$ and require truncation. – BruceET Nov 19 '17 at 15:48