[R] binomial simulation

[Moshe Olshansky]

Thank you - I wasn't aware of this function.
One can even use lchoose which allows really huge
arguments (more than 2^1000)!

--- "Lucke, Joseph F" <Joseph.F.Lucke at uth.tmc.edu>
wrote:

> C is an R function for setting contrasts in a
> factor.  Hence the funky
> error message.
> ?C
> 
> Use choose() for your C(N,k)
> ?choose
> 
> choose(200,2)
> 19900
> 
> choose(200,100)
>  9.054851e+58 
> 
> N=200; k=100; m=50; p=.6; q=.95
>
choose(N,k)*p^k*(1-p)^(N-k)*choose(k,m)*q^m*(1-q)^(k-m)
> 6.554505e-41
> 
> -----Original Message-----
> From: r-help-bounces at stat.math.ethz.ch
> [mailto:r-help-bounces at stat.math.ethz.ch] On Behalf
> Of Moshe Olshansky
> Sent: Wednesday, August 15, 2007 2:06 AM
> To: sigalit mangut-leiba; r-help
> Subject: Re: [R] binomial simulation
> 
> No wonder that you are getting overflow, since
> gamma(N+1) = n! and 200! > (200/e)^200 > 10^370.
> There exists another way to compute C(N,k). Let me
> know if you need this
> and I will explain to you how this can be done.
> But do you really need to compute the individual
> probabilities? May be
> you need something else and there is no need to
> compute the individual
> probabilities?
> 
> Regards,
> 
> Moshe.
> 
> --- sigalit mangut-leiba <smangut at gmail.com> wrote:
> 
> > Thank you,
> > I'm trying to run the joint probabilty:
> > 
> > C(N,k)*p^k*(1-p)^(N-k)*C(k,m)*q^m*(1-q)^(k-m)
> > 
> > and get the error: Error in C(N, k) : object not
> interpretable as a 
> > factor
> > 
> > so I tried the long way:
> > 
> > gamma(N+1)/(gamma(k+1)*(gamma(N-k)))
> > 
> > and the same with k, and got the error:
> > 
> > 1: value out of range in 'gammafn' in: gamma(N +
> 1)
> > 2: value out of range in 'gammafn' in: gamma(N -
> k) ....
> > 
> > Do you know why it's not working?
> > 
> > Thanks again,
> > 
> > Sigalit.
> > 
> > 
> > 
> > On 8/14/07, Moshe Olshansky
> <m_olshansky at yahoo.com>
> > wrote:
> > >
> > > As I understand this,
> > > P(T+ | D-)=1-P(T+ | D+)=0.05
> > > is the probability not to detect desease for a
> > person
> > > at ICU who has the desease. Correct?
> > >
> > > What I asked was whether it is possible to
> > mistakenly
> > > detect the desease for a person who does not
> have
> > it?
> > >
> > > Assuming that this is impossible the formula is
> > below:
> > >
> > > If there are N patients, each has a probability
> p
> > to
> > > have the desease (p=0.6 in your case) and q is
> the probability to 
> > > detect the desease for a person who
> > has
> > > it (q = 0.95 for ICU and q = 0.8 for a regular
> > unit),
> > > then
> > >
> > > P(k have the desease AND m are detected) = P(k
> have the desease)*P(m
> 
> > > are detected / k have
> > the
> > > desease) =
> > > C(N,k)*p^k*(1-p)^(N-k)*C(k,m)*q^m*(1-q)^(k-m)
> > > where C(a,b) is the Binomial coefficient "a
> above
> > b" -
> > > the number of ways to choose b items out of a
> > (when
> > > the order does not matter). You of course must
> > assume
> > > that N >= k >= m >= 0 (otherwise the probability
> > is
> > > 0).
> > >
> > > To generate such pairs (k infected and m
> detected)
> > you
> > > can do the following:
> > >
> > > k <- rbinom(N,1,p)
> > > m <- rbinom(k,1,q)
> > >
> > > Regards,
> > >
> > > Moshe.
> > >
> > > --- sigalit mangut-leiba <smangut at gmail.com>
> > wrote:
> > >
> > > > Hi,
> > > > The probability of false detection is: P(T+ |
> D-)=1-P(T+ | 
> > > > D+)=0.05.
> > > > and I want to find the joint probability
> > > > P(T+,D+)=P(T+|D+)*P(D+)
> > > > Thank you for your reply,
> > > > Sigalit.
> > > >
> > > >
> > > > On 8/13/07, Moshe Olshansky
> > <m_olshansky at yahoo.com>
> > > > wrote:
> > > > >
> > > > > Hi Sigalit,
> > > > >
> > > > > Do you want to find the probability P(T+ = t
> > AND
> > > > D+ =
> > > > > d) for all the combinations of t and d (for
> > ICU
> > > > and
> > > > > Reg.)?
> > > > > Is the probability of false detection (when
> > there
> > > > is
> > > > > no disease) always 0?
> > > > >
> > > > > Regards,
> > > > >
> > > > > Moshe.
> > > > >
> > > > > --- sigalit mangut-leiba <smangut at gmail.com>
> > > > wrote:
> > > > >
> > > > > > hello,
> > > > > > I asked about this simulation a few days
> > ago,
> > > > but
> > > > > > still i can't get what i
> > > > > > need.
> > > > > > I have 2 units: icu and regular. from icu
> I
> > want
> > > > to
> > > > > > take 200 observations
> > > > > > from binomial distribution, when
> probability
> > for
> > > > > > disease is: p=0.6.
> > > > > > from regular I want to take 300
> observation
> > with
> > > > the
> > > > > > same probability: p=0.6
> > > > > > .
> > > > > > the distribution to detect disease when
> > disease
> > > > > > occurred- *for someone from
> > > > > > icu* - is: p(T+ | D+)=0.95.
> > > > > > the distribution to detect disease when
> > disease
> > > > > > occurred- *for someone from
> > > > > > reg.unit* - is: p(T+ | D+)=0.8.
> > > > > > I want to compute the joint distribution
> for
> > > > each
> 
=== message truncated ===