[R] how to find "p" in binomial(n,p)

Charles C. Berry cberry at tajo.ucsd.edu
Thu Sep 20 05:17:28 CEST 2007


On Thu, 20 Sep 2007, Rolf Turner wrote:

>
> On 20/09/2007, at 11:04 AM, Mike Meredith wrote:
>
>>
>> I think the function you need is 'help.search'; try:
>>
>> help.search("binomial")
>>
>> and look for something obvious in the 'stats' package. A good deal
>> quicker
>> and easier than posting to an internet forum!
>>
>
> 	I don't think so.  I couldn't find anything useful under help.search
> ("binomial"), and
> 	I'm fairly comfortable with R.


After 20-odd years of R, Splus, and S, I'd say I am more enthusiastic than 
'comfortable', but trying 'help.search("binomial")' I get:

----
Help files with alias or concept or title matching 'binomial' using fuzzy 
matching:

[stuff omitted]

binom.test(stats)                   Exact Binomial Test


----

And after reading the help page:

> binom.test(15,20,conf=.95,alt='less') # Rolf's example

         Exact binomial test

data:  15 and 20
number of successes = 15, number of trials = 20, p-value = 0.994
alternative hypothesis: true probability of success is less than 0.5
95 percent confidence interval:
  0.0000000 0.8959192

...........^^^^^^^^^... matches Rolf's example

sample estimates:
probability of success
                   0.75

> binom.test(4,10,conf=.05,alt='le') # Ted's Example

         Exact binomial test

data:  4 and 10
number of successes = 4, number of trials = 10, p-value = 0.3770
alternative hypothesis: true probability of success is less than 0.5
5 percent confidence interval:
  0.0000000 0.2224411
...........^^^^^^^^^......... matches Ted's Example

sample estimates:
probability of success
                    0.4

So, it really can be found using 'help.search'.

Although in view of the 'Zen in the Art of Motorcycle Maintenance' effect 
(if it is only obvious if you know something, then it is not really 
obvious) I must note that one must recognize that "Pr(X <= k) <= alpha" is 
the confidence statement.

But of course, if one 'knew something' it is

 	help.search("beta")

that would be tried. :-)


HTH,

Chuck


>
> 	It seems to me that what is required here is to roll your own.
>
> 	I defined a function foo():
>
> 		foo <- function(p,n,k,alpha){pbinom(k,n,p)-alpha}
>
> 	Then I used uniroot to solve foo(p,n,k,alpha) = 0.
>
> 	E.g.
>
> 		xxx <- uniroot(foo,c(0,1),n=20,k=15,alpha=0.05)
>
> 	Then xxx$root is 0.895913.
>
> 	Check this out graphically:
>
> 		p <- seq(0,1,length=301)
> 		y <- foo(p,20,15,0.05)
> 		plot(p,y,type="l")
> 		abline(h=0)
> 		abline(v=xxx$root)
>
> 	It looks good.
>
> 				cheers,
>
> 					Rolf Turner
>> Cheers, Mike.
>>
>>
>> cathelf wrote:
>>>
>>> Dear all,
>>>
>>> I am trying a find the value "p" in binomial.
>>>
>>> X ~ Bin(n,p)
>>>
>>> I want to find the value "p", so that Pr(X <= k) <= alpha
>>>
>>> Here, n, k and alpha are known. n, k are integers. alpha is
>>> between (0,1).
>>>
>>> Thanks a lot!
>>>
>>> Catherine
>>>
>
>
> ######################################################################
> Attention:\ This e-mail message is privileged and confidenti...{{dropped}}
>
> ______________________________________________
> R-help at r-project.org 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.
>

Charles C. Berry                            (858) 534-2098
                                             Dept of Family/Preventive Medicine
E mailto:cberry at tajo.ucsd.edu	            UC San Diego
http://famprevmed.ucsd.edu/faculty/cberry/  La Jolla, San Diego 92093-0901



More information about the R-help mailing list