Tutorial 8b - Comparing two populations (Bayesian)

23 April 2011

> library(R2jags)
> library(ggplot2)
> library(grid)
> murray_opts <- opts(panel.grid.major = theme_blank(), 
+     panel.grid.minor = theme_blank(), 
+     panel.border = theme_blank(), panel.background = theme_blank(), 
+     axis.title.y = theme_text(size = 15, 
+         vjust = 0, angle = 90), axis.text.y = theme_text(size = 12), 
+     axis.title.x = theme_text(size = 15, 
+         vjust = -1), axis.text.x = theme_text(size = 12), 
+     axis.line = theme_segment(), plot.margin = unit(c(0.5, 
+         0.5, 1, 2), "lines"))

Bayesian t-test

Furness & Bryant (1996) studied the energy budgets of breeding northern fulmars (Fulmarus glacialis) in Shetland. As part of their study, they recorded the body mass and metabolic rate of eight male and six female fulmars.

Format of furness.csv data files
SEX METRATE BODYMASS MALE 2950 875 FEMALE 1956 765 MALE 2308 780 MALE 2135 790 MALE 1945 788 SEX Sex of breeding northern fulmars (Fulmarus glacialis) METRATE Metabolic rate (hJ/day) BODYMASS Body mass (g)

Open the furness data file.

> furness <- read.table("../downloads/data/furness.csv", 
+     header = T, sep = ",")
> head(furness)

           SEX METRATE BODYMASS
1 Male            2950      875
2 Female          1956      635
3 Male            2309      765
4 Male            2136      780
5 Male            1946      790
6 Female          1490      635

1. The researchers were interested in testing whether there is a difference in the metabolic rate of male and female breeding northern fulmars. In light of this, they are likely to want to answer the following:
   1. What is the difference (if any) in metabolic rate between males and females?
      
   2. What is the probability that males have a higher metabolic rate than females?
      

From a frequentist perspective, the appropriate statistical test for testing the null hypothesis that the means of two independent populations are equal is a t-test. We will fit a simple Bayesian linear model:

metrate_i = α + β_j(i)*sex_i + ε_i where ε ∼ N(0,σ²)