myBernBeta1 = function( priorShape , dataVec  ,mycex=1.4) {
# Bayesian updating for Bernoulli likelihood and beta prior.
# Input arguments:
#   priorShape
#     vector of parameter values for the prior beta distribution.
#   dataVec
#     vector of 1's and 0's.
#   credMass
#     the probability mass of the HDI.
# Output:
#   postShape
#     vector of parameter values for the posterior beta distribution.
# Graphics:
#   Creates a three-panel graph of prior, likelihood, and posterior
#   with highest posterior density interval.


# Check for errors in input arguments:
if ( length(priorShape) != 2 ) {
   stop("priorShape must have two components.") }
if ( any( priorShape <= 0 ) ) {
   stop("priorShape components must be positive.") }
if ( any( dataVec != 1 & dataVec != 0 ) ) {
   stop("dataVec must be a vector of 1s and 0s.") }

# Rename the prior shape parameters, for convenience:
a = priorShape[1]
b = priorShape[2]
# Create summary values of the data:
z = sum( dataVec == 1 ) # number of 1's in dataVec
N = length( dataVec )   # number of flips in dataVec
# Compute the posterior shape parameters:
postShape = c( a+z , b+N-z )

# Now plot everything:
# Construct grid of theta values, used for graphing.
binwidth = 0.005 # Arbitrary small value for comb on Theta.
Theta = seq( from = binwidth/2 , to = 1-(binwidth/2) , by = binwidth )
# Compute the prior at each value of theta.
pTheta = dbeta( Theta , a , b )
# Compute the likelihood of the data at each value of theta.
pDataGivenTheta = Theta^z * (1-Theta)^(N-z)
# Compute the posterior at each value of theta.
pThetaGivenData = dbeta( Theta , a+z , b+N-z )

# Open a window with three panels.
layout( matrix( c( 1,2,3 ) ,nrow=3 ,ncol=1 ,byrow=FALSE ) ) # 3x1 panels
par( mar=c(3,3,1,0) , mgp=c(2,1,0) , mai=c(0.5,0.5,0.3,0.1) ) # margin specs
maxY = max( c(pTheta,pThetaGivenData) ) # max y for plotting
# Plot the prior.
plot( Theta , pTheta , type="l" , lwd=3 ,
      xlim=c(0,1) , ylim=c(0,maxY) , cex.axis=1.2 ,
      xlab=bquote(theta) , ylab=bquote(p(theta)) , cex.lab=1.5 ,
      main="Prior" , cex.main=1.5 )
if ( a > b ) { textx = 0 ; textadj = c(0,1) } 
else { textx = 1 ; textadj = c(1,1) }
text( textx , 1.0*max(pThetaGivenData) ,
      bquote( "beta(" * theta * "|" * .(a) * "," * .(b) * ")"  ) ,
      cex=mycex ,adj=textadj )
# Plot the likelihood: p(data|theta)
plot( Theta , pDataGivenTheta , type="l" , lwd=3 ,
      xlim=c(0,1) , cex.axis=1.2 , xlab=bquote(theta) ,
      ylim=c(0,1.1*max(pDataGivenTheta)) ,
      ylab=bquote( "p(D|" * theta * ")" ) ,
      cex.lab=1.5 , main="Likelihood" , cex.main=1.5 )
if ( z > .5*N ) { textx = 0 ; textadj = c(0,1) }
else { textx = 1 ; textadj = c(1,1) }
text( textx , 1.0*max(pDataGivenTheta) , cex=mycex ,
      bquote( "Data: z=" * .(z) * ",N=" * .(N) ) ,adj=textadj )
# Plot the posterior.
plot( Theta , pThetaGivenData  ,type="l" , lwd=3 ,
      xlim=c(0,1) , ylim=c(0,maxY) , cex.axis=1.2 ,
      xlab=bquote(theta) , ylab=bquote( "p(" * theta * "|D)" ) ,
      cex.lab=1.5 , main="Posterior" , cex.main=1.5 )
if ( a+z > b+N-z ) { textx = 0 ; textadj = c(0,1) }
else { textx = 1 ; textadj = c(1,1) }

text( textx , 1.00*max(pThetaGivenData) , cex=mycex ,
      bquote( "beta(" * theta * "|" * .(a+z) * "," * .(b+N-z) * ")"  ) ,
      adj=textadj )

return( postShape )
} # end of function