MTH6102代写、代写R设计编程

日期：2023-01-06 01:05

Bayesian Statistical Methods

Queen Mary, Autumn 2022

1

Contents

1 Likelihood 4

1.1 Maximum likelihood estimation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5

1.2 Standard error . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8

2 Bayesian inference 10

2.1 Bayes＊ theorem . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10

2.2 Bayes＊ theorem and Bayesian inference . . . . . . . . . . . . . . . . . . . . . . . . . . 11

2.3 Conjugate prior distributions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13

2.4 Point estimates and credible intervals . . . . . . . . . . . . . . . . . . . . . . . . . . . 17

3 Specifying prior distributions 21

3.1 Informative prior distributions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21

3.2 Less informative prior distributions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22

4 Markov chain Monte Carlo methods 23

4.1 The Metropolis algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23

5 Predictive distributions 28

5.1 Simulating the posterior predictive distribution . . . . . . . . . . . . . . . . . . . . . 30

6 Hierarchical models 31

6.1 Inference for hierarchical models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32

6.2 Advantages of hierarchical models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35

7 Tests and model comparisions 36

2

Introduction

So far at Queen Mary, the statistics modules have taught the classical or frequentist approach which is

based on the idea that probability represents a long run limiting frequency. In the Bayesian approach,

any uncertain quantity is described by a probability distribution, and so probability represents a

degree of belief in an event which is conditional on the knowledge of the person concerned.

This course will introduce you to Bayesian statistics. These notes are self-contained but you may

want to read other accounts of Bayesian statistics as well. A useful introductory textbook is:

Bayesian Statistics: An Introduction (4th ed.) by P M Lee. (This is available as an e-book

from the library)

Parts of the following are also useful:

Bayesian Inference for Statistical Analysis by G E P Box and G C Tiao.

Bayesian Data Analysis by A Gelman, J B Carlin, H S Stern and D B Rubin.

Probability and Statistics from a Bayesian viewpoint (Vol 2) by D V Lindley.

3

Chapter 1

Likelihood

First we review the concept of likelihood, which is essential for Bayesian theory, but can also be used

in frequentist methods. Let y be the data that we observe, which is usually a vector. We assume

that y was generated by some probability model which we can specify. Suppose that this probability

model depends on one or more parameters 牟, which we want to estimate.

Definition 1.1. If the components of y are continuous, then the likelihood is defined as the joint

probability density function of y; if y is discrete, then the likelihood is defined as the joint probability

mass function of y. In either case, we denote the likelihood as

p(y | 牟).

Example 1.1. Let y = y1, . . . , yn be a random sample from a normal distribution with unknown

parameters 米 and 考. Then 牟 is the vector (米, 考), and the likelihood is the joint probability density

function

p(y | 米, 考) =

n﹊

i=1

?(yi | 米, 考).

Note that

﹊

here is the symbol for the product:

n﹊

i=1

ai = a1 ℅ a2 ℅ ﹞ ﹞ ﹞ ℅ an.

is the normal probability density function with parameters 米 and 考

(yi | 米, 考) = e(yi?米)2/2考2﹟2羽 考2.

4

Example 1.2. Suppose we observe k successes in n independent trials, where each trial has proba-

bility of success q. Now 牟 = q, the observed data is k, and the likelihood is the binomial probability

mass function

p(k | q) =(nk)

qk(1 q)nk.

It is also possible to construct likelihoods which combine probabilities and probability density func-

tions, for example if the observed data contains both discrete and continuous components. Alter-

natively, probabilities may appear in the likelihood if continuous data is only observed to lie within

some interval.

Example 1.3. Assume that the time until failure for a certain type of light bulb is exponentially

distributed with parameter 竹, and we observe n bulbs, with failure times t = t1, . . . , tn.

The likelihood contribution for a single observation ti is the exponential probability density function

竹e竹ti .

So the overall likelihood is the joint probability density function

p(t | 竹) =

n﹊

i=1

竹e?竹ti .

Suppose instead that we observe the failure time for the first m light bulbs with m < n, but for the

remaining n ?m bulbs we only observe that they have not failed by time ti. Then for i ≒ m, the

likelihood contributions are as before.

For i > m, the likelihood is the probability of what we have observed. Denoting the random variable

for the failure time by Ti, we have observed that Ti > ti, so the likelihood contribution is

p(Ti > ti) = e竹ti .

This is because the cumulative distribution function is p(Ti ≒ ti) = 1? e?竹ti .

Hence the overall likelihood is now

p(t | 竹) =

m﹊

i=1

竹e竹ti

n﹊

i=m+1

e竹ti .

1.1 Maximum likelihood estimation

In example 1.1, the parameters 米 and 考 of the normal distribution which generated the data are

known as population parameters. An estimator is defined as a function of the observed data which

5

we use as an estimate of a population parameter. For example the sample mean and variance may

be used as estimators of the population mean and variance, respectively.

To use the likelihood to estimate parameters 牟, we find the vector 牟? which maximizes the likelihood

p(y | 牟), given the data x that we have observed. This is the method of maximum likelihood, and

the estimator 牟? is known as the maximum likelihood estimator, or MLE.

Usually the parameters 牟 are continuous, and those are the only examples we cover, but the idea of

likelihood also makes sense if the unknown quantity is discrete.

When finding the MLE it is usually more convenient to work with the log of the likelihood: as the

log is a monotonically increasing function, the same 牟 will maximize p(y | 牟) and log(p(y | 牟)). The

log-likelihood is denoted by

?(牟; y) = log(p(y | 牟)).

Since the likelihood is typically a product of terms for independent observations, the log-likelihood

is a sum of terms, so using the log greatly simplifies finding the derivatives in order to find the

maximum.

Returning to the binomial example 1.2, the log-likelihood is

We also do not cover confidence intervals, but we do cover the Bayesian version, which are called

credible intervals.

9

Chapter 2

Bayesian inference

2.1 Bayes＊ theorem

Bayes＊ theorem is a formula from probability theory that is central to Bayesian inference. It is

named after Rev. Thomas Bayes, a nonconformist minister who lived in England in the first half of

the eighteenth century. The theorem states that:

Theorem 2.1. Let ? be a sample space and A1, A2, . . . , Am be mutually exclusive and exhaustive

events in ? (i.e. Ai ﹎ Aj = ?, i ?= j ,﹍k i=1Ai = ?; the Ai form a partition of ?.) Let B be any event

with p(B) > 0. Then

p(Ai | B) = p(Ai)p(B | Ai)

p(B)

=

p(Ai)p(B | Ai)﹉m

j=1 p(Aj)p(B | Aj)

The proof follows from the definition of conditional probabilities and the law of total probabilities.

Example 2.1. Suppose a test for an infection has 90% sensitivity and 95% specificity, and the

proportion of the population with the infection is q = 1/2000. Sensitivity is the probability of

detecting a genuine infection, and specificity is the probability of being correct about a non-infection.

So p(+ve test | infected) = 0.9 and p(-ve test | not infected) = 0.95. What is the probability that

someone who tests positive is infected?

Let the events be as follows:

B: test positive

A1: is infected

A2: is not infected

We want to find p(A1 | B). The probabilities we have are p(A1) = 1/2000, p(B | A1) = 0.9 and

p(B | A2) = 1? (BC | A2) = 1? 0.95 = 0.05.

10

Applying Bayes＊ theorem,

p(A1 | B) = p(A1)p(B | A1)﹉2

j=1 p(Aj)p(B | Aj)

p(A1)p(B | A1) = 1/2000℅ 0.9 = 0.00045

p(A2)p(B | A2) = (1? 1/2000)℅ 0.05 = .05

Hence

p(A1 | B) = 0.00045

0.00045 + .05

= 0.0089

So there is a less than 1% chance that the person is infected if they test positive.

Bayes＊ Theorem is also applicable to probability densities.

Theorem 2.2. Let X , Y be two continuous r.v.＊s (possibly multivariate) and let f(x, y) be the joint

probability density function (pdf), f(x | y) the conditional pdf etc. then

f(x | y) = f(y | x) f(x)

f(y)

=

f(y | x) f(x)÷

f(y | x∩) f(x∩) dx∩

Alternatively, Y may be discrete, in which case f(y | x) and f(y) are probability mass functions.

2.2 Bayes＊ theorem and Bayesian inference

In the Bayesian framework, all uncertainty is specified by probability distributions. This includes

uncertainty about the unknown parameters. So we need to start with a probability distribution for

the parameters p(牟).

We then update the probability distribution for 牟 using the observed data y. This updating is done

using Bayes＊ theorem

p(牟 | y) = p(牟) p(y | 牟)

p(y)

.

p(y | 牟) is the likelihood for parameters 牟 given the observed data y.

We don＊t normally need to find the normalizing constant p(y), which is given by

p(y) =

÷

p(牟) p(y | 牟) d牟 or

﹉

牟

p(牟) p(y | 牟)

11

So the procedure is as follows:

Start with a distribution p(牟) - this is known as the prior distribution.

Combine the prior distribution with the likelihood p(y | 牟) using Bayes＊ theorem.

The resulting probability distribution p(牟 | y) is known as the posterior distribution.

We base our inferences about 牟 on this posterior distribution.

The use of Bayes＊ theorem can be summarized as

Posterior distribution ≦ prior distribution ℅ likelihood

Example 2.2. Suppose a biased coin has probability of heads q, and we observe k heads in n

independent coin tosses. We saw the binomial likelihood for this problem:

p(k | q) =

(

n

k

)

qk(1? q)n?k

For Bayesian inference, we need to specify a prior distribution for q. As q is a continuous quantity

between 0 and 1, the family of beta distributions is a reasonable choice.

If X ‵ Beta(汐 , 汕), its probability density function is

f(x) =

x汐?1(1? x)汕?1

B(汐 , 汕)

, 0 ≒ x ≒ 1

where B is the beta function and 汐 , 汕 > 0 are parameters.

B(汐 , 汕) =

÷ 1

0

x汐?1(1? x)汕?1 dx, also B(汐 , 汕) = 忙(汐)忙(汕)

忙(汐 + 汕)

.

The uniform distribution on [0, 1] is a special case of the beta distribution with 汐 = 1, 汕 = 1.

Combining prior distribution and likelihood, we have the posterior distribution p(q | k) given by

p(q | k) ≦ p(q) p(k | q) = q

汐?1(1? q)汕?1

B(汐 , 汕)

(

n

k

)

qk(1? q)n?k ≦ qk+汐?1(1? q)n?k+汕?1

We can recognize this posterior distribution as being proportional to the pdf of a Beta(k+汐 , n?k+汕)

distribution. Hence we do not need to explicitly work out the normalizing constant for p(q | k), we

can immediately see that q | k has a Beta(k + 汐 , n? k + 汕) distribution.

12

For a Bayesian point estimate for q, we summarize the posterior distribution, for example by its

mean. A Beta(汐 , 汕) rv has expected value

汐

汐 + 汕

. Hence our posterior mean for q is

E(q | k) = k + 汐

n+ 汐 + 汕

.

Recall that the maximum likelihood estimator is

q? =

k

n

.

So for large values of k and n, the Bayesian estimate and MLE will be similar, whereas they differ

more for smaller sample sizes.

In the special case that the prior distribution is uniform on [0, 1], the posterior distribution is Beta(k+

1, n? k + 1). and the posterior mean value for q is

E(q | k) = k + 1

n+ 2

.

2.3 Conjugate prior distributions

In the binomial example 2.2 with a beta prior distribution, we saw that the posterior distribution is

also a beta distribution. When we have the same family of distributions for the prior and posterior

for one or more parameters, this is known as a conjugate family of distributions. We say that the

family of beta distributions is conjugate to the binomial likelihood.

The binomial likelihood is

p(牟1, 牟2, 牟3 | y) d牟2 d牟3

In practical Bayesian inference, Markov Chain Monte Carlo methods, covered later in these notes,

are the most common means of approximating the posterior distribution. These methods produce an

approximate sample from the joint posterior density, and once this is done the marginal distribution

of each parameter is immediately available.

Example 2.6. Generating samples from the posterior distribution may also be helpful even if we

can calculate the exact posterior distribution. Suppose that the data are the outcome of a clinical

trial of two treatments for a serious illness, the number of deaths after each treatment. Let the data

be ki deaths out of ni patients, i = 1, 2 for the two treatments, and the two unknown parameters

are q1 and q2, the probability of death with each treatment. Assuming a binomial model for each

outcome, and independent Beta(汐i, 汕i) priors for each parameter, the posterior distribution is

qi | ki ‵ Beta(ki + 汐i, ni ? ki + 汕i), i = 1, 2.

We have independent prior distributions and likelihood, so the posterior distributions are also inde-

pendent.

p(q1, q2 | k1, k2) = p(q1 | k1) p(q2 | k2) ≦ p(k1 | q1)p(q1) p(k2 | q2)p(q2)

However, it is useful to think in terms of the joint posterior density for q1 and q2, as then we can

make probability statements involving both parameters. In this case, one quantity of interest is the

probability P (q2 < q1), i.e. does the second treatment have a lower death rate than the first. To

find this probability, we need to integrate the joint posterior density over the relevant region, which

is not possible to do exactly when it is a product of beta pdfs.

We can approximate the probability by generating a sample of (q1, q2) pairs from the joint density.

To generate the sample, we just need to generate each parameter from its beta posterior distribution,

which can be done in R using the rbeta command. Then once we have the sample, we just count

what proportion of pairs has q2 < q1 to estimate P (q2 < q1).

20

Chapter 3

Specifying prior distributions

The posterior distribution depends on both the observed data via the likelihood, and also on the

prior distribution. So far, we have taken the prior distribution as given, but now we look at how to

specify a prior.

3.1 Informative prior distributions

An informative prior distribution is one in which the probability mass is concentrated in some subset

of the possible range for the parameter(s), and is usually based on some specific information. There

may be other data that is relevant, and we might want to use this information without including all

previous data in our current model. In that case, we can use summaries of the data to find a prior

distribution.

Example 3.1. In example 2.3, the data was the lifetimes of light bulbs, t = (t1, . . . , tn), assumed to

be exponentially distributed with parameter 竹 (the failure rate, reciprocal of the mean lifetime). The

gamma distribution provides a conjugate prior distribution for 竹. Suppose that we had information

from several other similar types of light bulbs, which had observed failure rates r1, . . . , rk. Let the

sample mean and variance of these rates be m and v.

A Gamma(汐 , 汕) distribution has mean