Wednesday, September 5, 2007

R computations for a proportion using a beta prior

Today we talked about using a beta prior to learn about a proportion. Inference about p is done by use of the beta posterior distribution and prediction about future samples is done by means of the predictive distribution.

Here are the R computations for the cell-phone example. I'll first illustrate inference for the proportion p, and then I'll illustrate the use of the special function pbetap (in the LearnBayes package) to compute the beta-binomial predictive distribution to learn about the number of successes in a future sample.

> library(LearnBayes)
>
> a=6.8; b=2.5 # parameters of beta prior
> n=24; y=9 # sample size and number of yes's in sample
>
> a1=a+y; b1=b+n-y # parameters of beta posterior
>
> # I'll illustrate different types of inferences
>
> # a point estimate is given by the posterior mean
> a1/(a1+b1)
[1] 0.4744745
>
> # or you could find the posterior median
> qbeta(.5,a1,b1)
[1] 0.4739574
>
> # a 90% interval estimate is found by use of the 5th and 95th quantiles
> # of the beta curve
> qbeta(c(.05,.95),a1,b1)
[1] 0.3348724 0.6158472
>
> # we illustrate prediction by use of the function pbetap
> # suppose we take a future sample of 20
> # how many will be driving when using a cell phone?
>
> m=20; ys=0:m
> pred.probs=pbetap(c(a1,b1),m,ys)
>
> # display the predictive probabilities
>
> cbind(ys,pred.probs)
ys pred.probs
[1,] 0 7.443708e-05
[2,] 1 6.444416e-04
[3,] 2 2.897264e-03
[4,] 3 8.968922e-03
[5,] 4 2.139155e-02
[6,] 5 4.170364e-02
[7,] 6 6.884411e-02
[8,] 7 9.841322e-02
[9,] 8 1.236003e-01
[10,] 9 1.376228e-01
[11,] 10 1.365218e-01
[12,] 11 1.208324e-01
[13,] 12 9.524434e-02
[14,] 13 6.650657e-02
[15,] 14 4.075296e-02
[16,] 15 2.159001e-02
[17,] 16 9.665296e-03
[18,] 17 3.527764e-03
[19,] 18 9.889799e-04
[20,] 19 1.901993e-04
[21,] 20 1.891124e-05
>
> # what is the probability that there are at least
> # 10 cell phone drivers in my sample?
>
> sum(pred.probs*(ys>=10))
[1] 0.4958392
>

1 comment:

Seb said...

Dear Dr. Albert, in your book, the introductory example for the same model (for the proportion of sleepers), the parameters for the prior are said to be obtained by trial and error. I wonder if there is an analytical way to do it, using the properties of the beta distribution. Thanks.