Prior distribution to the binomial distribution probability distributions urn model I have an infinite population with unknown mean of successes and failures. I'm drawing 400 times from the population and get 400 successes. Now I want to generate random estimates for the true mean of the population from which I have drawn the 400 successes. Can anyone tell me which probability distribution I have to use, respectively how the function to use in R would look like?
I thought that
rbeta(1,400+1,400-400+1)

might be the right function, but even if I perform this one 10 000 000 times I never get the result 1 (and the same accounts for rbeta(10000000,0+1,400-0+1) and the result 0). So I ask myself why isn't it possible to have a population which consists only of successes or only of failures?   
 A: I am under the impression your spirit is a little confused. You have to distinguish three things:


*

*the real experiment

*the mathematical (probabilistic/statistical) model for the real experiment

*the computational aspects of the statistical inference with the mathematical model


For the first point, you claim your population is infinite. So 1) there's no real experiment and your question is purely related to the mathematical model; or 2) you have a real experiment with a huge population and the consideration of an infinite population actually is an assumption of the model (a model can only approximate the real experiment).
You assume a Beta posterior distribution for the proportion of successes. This distribution is a mathematical object, it is a part of the mathematical modeling of your problem. Firstly, a Beta random variable never takes the value $0$ - this is a matter of probability theory; more rigorously we should say that it almost never equals $0$, but the "almost never" probabilistic notion aims to modelize the "never" notion in the real world. Secondly, you use a computer to generate a Beta variable: the numbers you generate are not random, but pseudo-random, and everything is done by the programmer in order that $0$ and $1$ will never be generated (*).
Finally, I suggest you take a look at the probability to sample a success in your infinite population. To do so, you could use the Beta-binomial distribution available in the VGAM package or evaluate this probability using simulations:
> n <- 1000000
> sims <- rep(NA,n)
> for(i in 1:n){
+ p <- rbeta(1,401,1)
+ sims[i] <- rbinom(1,1,p)
+ }
> mean(x)
[1] 0.99754

It is not $1$, why ? Assume you sample $2$ individuals and you get $2$ successes in a population made of $10^6$ individuals. Would you infer that success occurs for every individual of the population ? So what about $400$ successes in an infinite population ?
(*) This point is wrong. See comments below. But for instance a computer never generates the values $0$ and $1$ for a uniform distribution on $[0,1]$.
