if (FALSE) #Make Romeo slide (Bayesian update, continuous-continuous 1)
{
  makepdf = FALSE
  theta1 = seq(0,.4,.01)
  theta2 = seq(.4,1,.01)
  pr1 = 0*theta1 + 1
  pr2 = 0*theta2 + 1
  T = 21/8
  post1 = 0*theta1
  post2 = 1/(T*theta2^3)
  if (makepdf)
  {
    fname = '../pdftex/img/class21-1.pdf'
    pdf(file=fname, width=3, height=2)
    par(mar=c(2,2,1.7,.3))
  } 
  title = expression(paste('Prior and posterior for ',theta))
  plot(c(0,1),c(0,max(post2)), type='n', xlab='', ylab='',
       main=title,cex.main=1.4,cex.lab=1.3, cex.axis = 1)
  lines(c(theta1,theta2), c(pr1,pr2), col='red', lwd=2)
  lines(c(theta1,theta2), c(post1,post2), col='cyan', lwd=2)
  abline(h=0)
  if (makepdf)
    dev.off()
}

if (FALSE) #Make waiting times 1 slide (Bayesian update, continuous-continuous 2)
{
  makepdf = FALSE
  lambda = seq(0,2.5,.01)
  pr = 2*exp(-2*lambda)
  T = 4/125
  post = (2/T)*lambda^2*exp(-5*lambda)
  if (makepdf)
  {
    fname = '../pdftex/img/class21-2.pdf'
    pdf(file=fname, width=3, height=2)
    par(mar=c(2,2,1.7,.3))
  } 
  title = expression(paste('Prior and posterior for ',lambda))
  plot(c(0,max(lambda)),c(0,max(pr)), type='n', xlab='', ylab='',
       main=title,cex.main=1.4,cex.lab=1.3, cex.axis = 1)
  lines(lambda, pr, col='red', lwd=2)
  lines(lambda, post, col='cyan', lwd=2)
  abline(h=0)
  if (makepdf)
    dev.off()
}

if (FALSE) #ANOVA question
{
  #DATA ----
  T1 = c(6,8,4,5,3,4)
  T2 = c(8,12,9,11,6,8)
  T3 = c(13,9,11,8,7,12)
  
  #One way ANOVA using aov() ----
  treat = c(rep('T1',length(T1)),rep('T2',length(T2)),rep('T3',length(T3)))
  recovery = c(T1,T2,T3)
  data.recovery = data.frame(treat,recovery)
  print(data.recovery)
  aov.data = aov(recovery~treat,data=data.recovery)  #do the analysis of variance
  summary(aov.data)                          #show the summary table
  print(model.tables(aov.data,"means"),digits=3)  #report the means and the number of subjects/cell
  boxplot(recovery~treat,data=data.recovery)#graphical summarydatafilename
  
  # One way ANOVA by hand ---- 
  # Make sure all the groups are the same size 
  if (length(T1) != length(T2) || length(T2) != length(T3))
    stop("lengths are not equal", call. = F)
  m = length(T1)
  n = 3
  group.means = c(mean(T1),mean(T2), mean(T3))
  grandmean = mean(group.means)
  group.vars = c(var(T1),var(T2), var(T3))
  betweenGroupVar = m*var(group.means)
  aveWithinGroupVar = mean(group.vars)
  fstat = betweenGroupVar/aveWithinGroupVar
  p = 1-pf(fstat, n-1,n*(m-1))
  print(paste('Group means',toString(group.means)))
  print(paste('Grandmean',grandmean))
  print(paste('Group vars', toString(group.vars)))
  print(paste('Between group var', betweenGroupVar))
  print(paste('Within group var', aveWithinGroupVar))
  print(paste('fstat',fstat))
  print(paste('p',p))
  
}