MATH3871/MATH5960
Bayesian Inference and Computation
Lab 1 Exercises
These exercises provide some practice in performing basic Bayesian analyses. Some involve
doing some algebra, others involve working with R. There is no requirement to do the exercises
in order – just attempt the questions that interest you the most. Outline solutions are available
in a separate file.
1) Water consumption
Water consumption
Water use
More than 5 million ML of water is used in Victoria each year, 90% from surface
water and 10% from groundwater. The majority of Victoria’s water resource is
used for irrigation (78%), while urban uses (both metropolitan and regional)
account for 17% of Victoria’s water consumption.
Efficiency of water use2
Making the best use of our water resource
CSIRO has provided a perspective on water use which is related to the
contribution it makes to Gross National Expenditure (GNE). The prices we pay
for agricultural products tend not to take into account the full environmental
costs of production. When water inputs are considered there is variation in
the amount of water used to produce various commodities. For example,
rice production uses more than 8,000 litres of water for every dollar of GNE,
compared to cereal production which uses around 600 litres to produce the
same value of product.2
Tr nds in water consumption
While time-series water consumption data are not available for regional Victoria,
records for Melbourne show that per capita water consumption grew steadily
between the 1940s and 1960s, with strong increases in the 1970s. Since the
1980s water consumption has been influenced by drought, and associated water
restrictions, as well as by conservation and efficiency incentives and market
reforms designed to reduce water consumption over the longer term.3
Average daily per capita water use4
Melbourne 1940-2004*
* NOTE: Figure for 2003-04 is forecasted estimation
Consumptive uses of water in Victoria1
2003-04
BTRc^a ;XcaTb^UfPcTa_Ta^U6a^bb=PcX^]P[4g_T]SXcdaT
AXRT
2^cc^]
3PXahBdVPaRP]T
3PXah_a^SdRcb
1TTU
5[^daRTaTP[U^^Sb
BTaeXRTbc^PVaXRd[cdaT
'#'
% #
#$!
!#&
'$!
&"
% '
$'#
Sources 1DSE 2005 State Water Report 2CSIRO and University of Sydney 2005 Balancing Act. A Triple Bottom Line Analysis of the Australian Economy 3Victorian Government 2004 Securing Our Water Future
Together 4Melbourne Water 2005 A Dry History
&

8aaXVPcX^]
AdaP[3^\TbcXRP]SBc^RZ
ATVX^]P[DaQP]
&'
$
1 9
4 0
1 9
4 2
1 9
4 4
1 9
4 6
1 9
4 8
1 9
5 0
1 9
5 2
1 9
5 4
1 9
5 6
1 9
5 8
1 9
6 0
1 9
6 2
1 9
6 4
1 9
6 6
1 9
6 8
1 9
7 0
1 9
7 2
1 9
7 4
1 9
7 6
1 9
7 8
1 9
8 0
1 9
8 2
1 9
8 4
1 9
8 6
1 9
8 8
1 9
9 0
1 9
9 2
1 9
9 4
1 9
9 6
1 9
9 8
2 0
0 0
2 0
0 2
2 0
0 4
Financial Year Ending
L i
t r
e s
p
e r
p
e r
s o
n
p e
r
d a
y
550
500
450
400
350
300
250
200
150
100
50
0
Drought Years
Water Restrictions
Mid 1990s:
COAG Water
Reforms - Pricing
measures introduced
86
NOTE: These data are for
Australia overall. Cotton is
not grown in Victoria and rice
represents less than 1% of
the State’s grain production
(as at 2001 Agricultural
Census)
In the Melbourne average daily per capita water use analyis, we modelled the discrete observa-
tions x1, . . . , xn as independent draws from a Poisson(θ) distribution. Assuming a Gamma(α, β)
prior, which has a density function of
pi(θ) =
βα
Γ(α)
θα−1 exp(−βθ), for α, β > 0,
we computed the posterior as a Gamma (α +
∑n
i=1 xi, β + n) distribution.
(a) Given that n = 65,
∑
i xi = 24, 890 and with prior parameters α = 1, β = 0.01, compute
a point estimate (i.e. the posterior mean) and a 95% central credible interval for θ. Note,
you will need to compute the credible interval numerically in R (hint: use the R command
qgamma).
(b) Draw a sample of sizeN = 500 directly from the posterior distribution (see the R command
rgamma), and obtain Monte Carlo estimates of the lower and upper values of the 95%
credible interval for θ.
1
Repeat this 100 times, and produce a histogram for the distribution of each interval end-
point. Superimpose a point corresponding to the true interval endpoints. How accurate
is the Monte Carlo estimate? (R commands: hist, points).
Produce another pair of histograms, but this time use N = 2500 samples. How is the
precision of the Monte Carlo estimates affected?
(c) Draw samples directly from the posterior distribution. Use these to obtain samples from
the posterior predictive distribution, and plot this via a histogram. Superimpose the den-
sity of the algebraically computed negative binomial predictive distribution (R command:
dnbinom). Do the distributions coincide?
(d) What are the advantages/disadvantages of performing statistical analyses using the alge-
braically exact approach, and the Monte Carlo approximations?
2
2) Buffon’s Needle
One of the most famous simulation experiments is Buffon’s Needle, designed to calculate (not
very efficiently!) an estimate of pi. Imagine a grid of parallel lines with spacing d, on which a
needle of lenght ` ≤ d is dropped. We repeat this experiment n times and count the proportion
of times pˆ that the needle intersects with a line.
The rationale behind this is that if x is the distance from the centre of the needle to the
leftmost line, and if θ is the angle from the vertical, then under the assumption of random
needle throwing, we would have x ∼ U(0, d) and θ ∼ U(0, pi). Hence
p = Pr(needle intersects line)
=
1
pi
∫ pi
0
Pr(needle intersects|θ = φ)dφ
1
pi
∫ pi
0
(
2
d
× `
2
sinφ
)
dφ
=
2`
pid
.
Hence, an estimate of pi is pˆi = 2`
pd
a) Produce some code to simulate the Buffon’s Needle experiment, given the lengths ` and
d, and produce an estimate of pi. Plot the estimate of pi as the number of simulations, n,
increases.
b) A natural question is how to optimise the relative sizes of ` and d. Consider the variability
of 1
pˆi
.
Now npˆ ∼ Bin(n, p), so Var(pˆ) = p(1− p)/n. Show that Var(1/pˆi) = 1
npi2
(
pi
2ρ
− 1
)
where
ρ = `/d. When is this minimised (for 0 ≤ ρ ≤ 1)?
c) By computing the estimate of pi 1000 times and computing the standard deviation, for
a range of values ρ = `/d, empirically demonstrate that your optimal value of ρ leads to
the smallest variability for pˆi.
There are a number of things which may (or may not!) improve the efficiency of this
experiment, including:
• using a grid of rectangles or squares;
3
• using a cross or other shape instead of a needle;
• using a needle of length greater than the grid separation.
The point is: simulation can be used to answer many interesting problems, but careful
design may be needed to achieve even moderate efficiency.
4