Learn R Programming

Bolstad (version 0.1-8)

normgcp: Bayesian inference on a normal mean with a general continuous prior

Description

Evaluates and plots the posterior density for mu, the mean of a normal distribution, with a general continuous prior on mu

Usage

normgcp(x, sigma.x, density = "uniform" , params = NULL, n.mu = 50, mu = NULL, mu.prior = NULL, ret = FALSE)

Arguments

x
a vector of observations from a normal distribution with unknown mean and known std. deviation.
sigma.x
the population std. deviation of the normal distribution
density
distributional form of the prior density can be one of: "normal", "unform", or "user".
params
if density = "normal" then params must contain at least a mean and possible a std. deviation. If a std. deviation is not specified then sigma.x will be used as the std. deviation of the prior. If density = "uniform" then params must contain a minimum and
n.mu
the number of possible mu values in the prior
mu
a vector of possibilities for the probability of success in a single trial. Must be set if density="user"
mu.prior
the associated prior probability mass. Must be set if density="user"
ret
if true then the likelihood and posterior are returned as a list.

Value

  • If ret is true, then a list will be returned with the following components:
  • likelihoodthe scaled likelihood function of x given mu and sigma.x
  • posteriorthe posterior probability of mu given x and sigma.x
  • muthe vector of possible mu values used in the prior
  • mu.priorthe associated probability mass for the values in mu

See Also

normdp normnp

Examples

Run this code
## generate a sample of 20 observations from a N(-0.5,1) population
x<-rnorm(20,-0.5,1)

## find the posterior density with a uniform U[-3,3] prior on mu
normgcp(x,1,params=c(-3,3))

## find the posterior density with a non-uniform prior on mu
mu<-seq(-3,3,by=0.1)
mu.prior<-rep(0,length(mu))
mu.prior[mu<=0]<-1/3+mu[mu<=0]/9
mu.prior[mu>0]<-1/3-mu[mu>0]/9
normgcp(x,1,density="user",mu=mu,mu.prior=mu.prior)

## find the CDF for the previous example and plot it
results<-normgcp(x,1,density="user",mu=mu,mu.prior=mu.prior,ret=TRUE)
cdf<-sintegral(mu,results$posterior,n.pts=length(mu),ret=TRUE)
plot(cdf,type="l",xlab=expression(mu[0])
             ,ylab=expression(Pr(mu<=mu[0])))

## use the CDF for the previous example to find a 95\%
## credible interval for mu. Thanks to John Wilkinson for this simplified code

lcb<-cdf$x[with(cdf,which.max(x[y<=0.025]))]
ucb<-cdf$x[with(cdf,which.max(x[y<=0.975]))]
cat(paste("Approximate 95% credible interval : ["
           ,round(lcb,4),"",round(ucb,4),"]
",sep=""))

## use the CDF from the previous example to find the posterior mean
## and std. deviation
dens<-mu*results$posterior 
post.mean<-sintegral(mu,dens)

dens<-(mu-post.mean)^2*results$posterior
post.var<-sintegral(mu,dens)
post.sd<-sqrt(post.var)

## use the mean and std. deviation from the previous example to find
## an approximate 95\% credible interval
lb<-post.mean-qnorm(0.975)*post.sd
ub<-post.mean+qnorm(0.975)*post.sd


cat(paste("Approximate 95% credible interval : ["
   ,round(lb,4),"",round(ub,4),"]
",sep=""))

Run the code above in your browser using DataLab