\chapter{Single parameter problems} This chapter reviews some of the key results for Bayesian inference of single parameter problems studied in MAS2903. \section{Prior and posterior distributions} Suppose we have data $\xvece$ which we model using the probability (density) function $f(\xvec|\theta)$, which depends on a single parameter~$\theta$. Once we have observed the data, $f(\xvec|\theta)$ is the \textit{likelihood function} for $\theta$ and is a function of~$\theta$ (for fixed $\xvec$) rather than of~$\xvec$ (for fixed~$\theta$). Also, suppose we have prior beliefs about likely values of $\theta$ expressed by a probability (density) function $\pi(\theta)$. We can combine both pieces of information using the following version of Bayes Theorem. The resulting distribution for $\theta$ is called the posterior distribution for~$\theta$ as it expresses our beliefs about $\theta$ {\it after} seeing the data. It summarises all our current knowledge about the parameter $\theta$. Using Bayes Theorem, the posterior probability (density) function for $\theta$ is \begin{align*} \pi(\theta|\xvec)&=\frac{\pi(\theta)\,f(\xvec|\theta)}{f(\xvec)} \\ \intertext{where} f(\xvec)&= \begin{cases} \int_\Theta\pi(\theta)\,f(\xvec|\theta)\,d\theta & \text{if $\theta$ is continuous}, \\ \\ \sum_\Theta\pi(\theta)\,f(\xvec|\theta) & \text{if $\theta$ is discrete}. \end{cases} \intertext{Also, as $f(\xvec)$ is not a function of $\theta$, Bayes Theorem can be rewritten as} \pi(\theta|\xvec)&\propto \pi(\theta)\times f(\xvec|\theta) \\ i.e.~~\text{Posterior}&\propto\text{Prior}\times\text{Likelihood}. \end{align*} \begin{example}~\\ Table~\ref{tab:botulism} shows some data on the number of cases of \label{ex:botulism} foodbourne botulism in England and Wales. It is believed that cases occur at random at a constant rate~$\theta$ in time (a Poisson process) and so can be modelled as a random sample from a Poisson distribution with mean~$\theta$. \begin{table}[ht] \begin{center} \begin{tabular}{|lcccccccc|} \hline Year & 1998 & 1999 & 2000 & 2001 & 2002 & 2003 & 2004 & 2005 \\ Cases & 2 & 0 & 0 & 0 & 1 & 0 & 2 & 1\\ \hline \end{tabular} \caption{Number of cases of foodbourne botulism in England and Wales, 1998--2005} \label{tab:botulism} \end{center} \end{table} \framebox[17cm]{ \parbox[t]{16cm}{ \vspace{0.4cm} An expert in the epidemiology of similar diseases was consulted to construct a suitable prior distribution for the rate ~$\theta$. \begin{itemize} \item [(a)] What might be an appropriate probability model for the rate $\theta$? Why is this termed the \textit{prior distribution}? \huge \ding{46} \normalsize \textcolor{blue}{A gamma distribution might be appropriate here as $\theta$ is a rate and must be positive; the gamma distribution has support over the positive real line. This is termed the prior distribution as it represents our beliefs about $\theta$ before observing any data.} %\vspace{3cm} \item [(b)] Suppose the epidemiologist tells you that she expects around two cases of foodbourne botulism per year in England and Wales. After further discussions with the epidemiologist, you are also able to elicit that $SD(\theta) \approx \sqrt{2}$. Using your answer to part (a) find your full prior distribution for $\theta$, and write down its density. \huge \ding{46} \normalsize \textcolor{blue}{Suppose our prior is $\theta \sim Ga(g, h)$. Then from the information given, we can use the formulae for the gamma mean and standard deviation to write down: \[ E(\theta) = \frac{g}{h} = 2 \quad \text{(using the expert's most likelt value for $\theta$)} \] and \[ SD(\theta) = \frac{\sqrt{g}}{h} = \sqrt{2}. \] Solving simulataneously for $g$ and $h$ gives $g=2$ and $h=1$, and so we have $\theta \sim Ga(2,1)$. Therefore: \[ \pi(\theta) = \frac{1^2 \theta^{2-1} e^{-1\theta}}{\Gamma(2)} = \propto \theta e^{-\theta}, \quad \theta>0. \] } %\vspace{7cm} \end{itemize} }} \framebox[17cm]{ \parbox[t]{16cm}{ \vspace{0.4cm} \begin{itemize} \item [(c)] Determine the posterior distribution for ~$\theta$. \huge \ding{46} \normalsize \textcolor{blue}{ We know that \[ \text{Posterior} \propto \text{Prior} \times \text{Likelihood}. \] We know from the last question that the prior density is $\pi(\theta) \propto \theta e^{-\theta}$. What about the Poisson likelihood? From the data provided, we can write down that \begin{eqnarray*} f(\boldsymbol{x}|\theta) &= &\frac{e^{-\theta}\theta^2}{2!} \times \frac{e^{-\theta}\theta^0}{0!} \times \cdot \times \frac{e^{-\theta}\theta^1}{1!} \\ &\propto e^{-8\theta} \theta^{6}. \end{eqnarray*} Therefore \begin{eqnarray*} \pi(\theta|\boldsymbol{x}) &\propto& \pi(\theta) \times f(\boldsymbol{x}|\theta) \\ &\propto& \theta e^{-\theta} \times \theta^{6} e^{-8\theta} \\ &\propto& \theta^{7}e^{-9\theta} \quad \theta>0 \end{eqnarray*} and we recognise this as being proportional to the density of a $Ga(8, 9)$ distribution; that is, $\theta|\boldsymbol{x} \sim Ga(8, 9)$.} \vspace{5cm} %\vspace{15cm} \end{itemize} }} \end{example} Thus the data have updated our beliefs about $\theta$ from a $Ga(2,1)$ distribution to a $Ga(8,9)$ distribution. Plots of these distributions are given in Figure~\ref{fig:botulismpost}, and Table~\ref{tab:botulismsum} gives a summary of the main changes induced by incorporating the data --- a $Ga(g,h)$ distribution has mean $g/h$, variance $g/h^2$ and mode $(g-1)/h$. \begin{figure}[ht] \begin{center} \includegraphics[width=0.8\textwidth]{botulism} \caption{Prior (dashed) and posterior (solid) densities for $\theta$} \label{fig:botulismpost} \end{center} \end{figure} Notice that, as the mode of the likelihood function is close to that of the prior distribution, the information in the data is consistent with that in the prior distribution. Also there is a reduction in variability from the prior to the posterior distributions. The similarity between the prior beliefs and the data has reduced the uncertainty we have about the rate~$\theta$ at which cases occur. \begin{table}[ht] \begin{center} \begin{tabular}{|l|c|c|c|} \cline{2-4} \multicolumn{1}{c|}{~}& Prior & Likelihood & Posterior \\ \hline $Mode(\theta)$ & 1.00 & 0.75 & 0.78 \\ $E(\theta)$ & 2.00 & -- & 0.89 \\ $SD(\theta)$ & 1.41 & -- & 0.31 \\ \hline \end{tabular} \caption{Changes in beliefs about $\theta$} \label{tab:botulismsum} \end{center} \end{table} \begin{example}~\\ Consider now the general case of Example~\ref{ex:botulism}: \label{ex:poisgam} suppose $X_i|\theta\sim Po(\theta)$, $i=1,2,\ldots,n$ (independent) and our prior beliefs about $\theta$ are summarised by a $Ga(g,h)$ distribution (with $g$ and $h$ known), with density \begin{equation} \label{eq:poprior} \pi(\theta)=\frac{h^g\,\theta^{g-1}e^{-h\theta}}{\Gamma(g)}, \quad\theta>0. \end{equation} Determine the posterior distribution for $\theta$. \framebox[17cm]{ \parbox[t]{16cm}{ \vspace{0.4cm} \huge \ding{46} \normalsize \vspace{3cm}}} \newpage \framebox[17cm]{ \parbox[t]{16cm}{ \vspace{0.4cm} \huge \ding{46} \normalsize \textit{Solution to Example 1.2 continued...} \textcolor{blue}{ The likelihood function for~$\theta$ is \begin{align} \label{eq:polik} f(\xvec|\theta)&=\prod_{i=1}^n \frac{\theta^{x_i} e^{-\theta}}{x_i!}, \quad\quad\theta>0 \notag \\ &\propto\theta^{n\bar x} e^{-n\theta}, \quad\quad\theta>0. \\ \intertext{Using Bayes Theorem, the posterior density function is} \pi(\theta|\xvec)&\propto\pi(\theta)\,f(\xvec|\theta) \notag \\ &\propto\frac{h^g\,\theta^{g-1}e^{-h\theta}}{\Gamma(g)}\times \theta^{n\bar x} e^{-n\theta},\quad\quad\theta>0 \notag \\ \text{i.e.}\quad\pi(\theta|\xvec) \label{eq:popost} &=k\theta^{g+n\bar x-1}e^{-(h+n)\theta},\quad\quad\theta>0 \end{align} where $k$ is a constant that does not depend on $\theta$. Therefore, the posterior density takes the form $k\theta^{G-1}e^{-H\theta}$, $\theta>0$ and so the posterior must be a gamma distribution. Thus we have $\theta|\xvec\sim Ga(G=g+n\bar x,H=h+n)$.} }} Summary: If we have a random sample from a $Po(\theta)$ distribution and our prior beliefs about $\theta$ follow a $Ga(g,h)$ distribution then, after incorporating the data, our (posterior) beliefs about $\theta$ follow a $Ga(g+n\bar x,h+n)$ distribution. The changes in our beliefs about $\theta$ are summarised in Table~\ref{tab:gam}, taking $g\geq 1$. \begin{table}[ht] \bigskip \begin{center} \begin{tabular}{|l|c|c|c|} \cline{2-4} \multicolumn{1}{c|}{~}& Prior & Likelihood & Posterior \\ \multicolumn{1}{c|}{~}& $\eqref{eq:poprior}$ & $\eqref{eq:polik}$ & $\eqref{eq:popost}$ \\ \hline $Mode(\theta)$ & $(g-1)/h$ & $\bar x$ & $(g+n\bar x-1)/(h+n)$ \\ $E(\theta)$ & $g/h$ & -- & $(g+n\bar x)/(h+n)$ \\ $SD(\theta)$ & $\sqrt{g}/h$ & -- & $\sqrt{g+n\bar x}/(h+n)$ \\ \hline \end{tabular} \caption{Changes in beliefs about $\theta$} \label{tab:gam} \end{center} \end{table} Notice that the posterior mean is greater than the prior mean if and only if the likelihood mode is greater than the prior mean, that is, \begin{equation*} E(\theta|\xvec)>E(\theta)\quad\iff\quad Mode_\theta\{f(\xvec|\theta)\}>E(\theta). \end{equation*} The standard deviation of the posterior distribution is smaller than that of the prior distribution if and only if the sample mean is not too large, that is \[ SD(\theta|\xvec)E(\mu)\quad\iff\quad Mode_\mu\{f(\xvec|\mu)\}>E(\mu). \end{equation*} Also, the standard deviation of the posterior distribution is smaller than that of the prior distribution. \end{example} \begin{example}~\\ The 18th century physicist Henry Cavendish made 23 experimental \label{ex:cavendish} determinations of the earth's density, and these data (in $g/cm^3$) are given below. \begin{center} \begin{tabular}{|cccccccc|} \hline 5.36&5.29&5.58&5.65&5.57&5.53&5.62&5.29\\ 5.44&5.34&5.79&5.10&5.27&5.39&5.42&5.47\\ 5.63&5.34&5.46&5.30&5.78&5.68&5.85&\\ \hline \end{tabular} \end{center} Suppose that Cavendish asserts that the error standard deviation of these measurements is $0.2\,g/cm^3$, and assume that they are normally distributed with mean equal to the true earth density~$\mu$. Using a normal prior distribution for $\mu$ with mean $5.41\,g/cm^3$ and standard deviation $0.4\,g/cm^3$, derive the posterior distribution for $\mu$. \framebox[17cm]{ \parbox[t]{16cm}{ \vspace{0.4cm} \huge \ding{46} \normalsize \textcolor{blue}{ From the data we calculate $\bar x=5.4848$ and $s=0.1882$. Therefore, the assumed standard deviation $\sigma=0.2$ is probably okay. We also have $\tau=1/0.2^2$, $b=5.41$, $d=1/0.4^2$ and $n=23$.} \vspace{2.5cm}}} \newpage \framebox[17cm]{ \parbox[t]{16cm}{ \vspace{0.4cm} \huge \ding{46} \normalsize \textit{Solution to Example 1.4 continued...} \textcolor{blue}{ Therefore, using Example~\ref{ex:normal}, the posterior distribution is $\mu|\xvec\sim N(B,1/D)$, where \begin{align*} B&=\frac{db+n\tau\bar x}{d+n\tau} =\frac{5.41/0.4^2+23\times 5.4848/0.2^2}{1/0.4^2+23/0.2^2} =5.4840 \intertext{and} D&=d+n\tau=\frac{1}{0.4^2}+\frac{23}{0.2^2}=\frac{1}{0.0415^2}. \end{align*} Therefore the posterior distribution is $\mu|\xvec\sim N(5.484,0.0415^2)$ and is shown in Figure~\ref{fig:cavendishplot}. } }} %%Hack for vspace %\vspace{-3.5cm} The actual mean density of the earth is $5.515\,g/cm^3$ (Wikipedia). We can determine the (posterior) probability that the mean density is within 0.1 of this value as follows. The posterior distribution is $\mu|\xvec\sim N(5.484,0.0415^2)$ and so \[ Pr(5.415<\mu<5.615|\xvec)=0.9510, \] calculated using the \verb$R$ command \verb$pnorm(5.615,5.484,0.0415)-pnorm(5.415,5.484,0.0415)$. Without the data, the only basis for determining the earth's density is via the prior distribution. Here the prior distribution is $\mu\sim N(5.4,0.4^2)$ and so the (prior) probability that the mean density is within 0.2 of the (now known) true value is \[ Pr(5.315<\mu<5.715)=0.1896, \] calculated using the \verb$R$ command \verb$pnorm(5.615,5.4,0.4)-pnorm(5.415,5.4,0.4)$. \begin{figure}[ht] \begin{center} \includegraphics[width=0.73\textwidth]{cavendish} \caption{Prior (dashed) and posterior (solid) densities for the earth's density} \label{fig:cavendishplot} \end{center} \end{figure} \end{example} \section{Different levels of prior knowledge} \subsection{Substantial prior knowledge} We have {\it substantial prior information} for $\theta$ when the prior distribution {\it dominates} the posterior distribution, that is $\pi(\theta|\xvec)\sim\pi(\theta)$. When we have substantial prior information there can be some difficulties: \begin{enumerate} \item the intractability of the mathematics in deriving the posterior distribution --- though with modern computing facilities this is less of a problem, \item the practical formulation of the prior distribution --- coherently specifying prior beliefs in the form of a probability distribution is far from straightforward. \end{enumerate} \subsection{Limited prior knowledge} When prior information about $\theta$ is limited, the pragmatic approach is to choose a distribution which makes the Bayes updating from prior to posterior mathematically straightforward, and use what prior information is available to determine the parameters of this distribution. For example \begin{itemize} \item Poisson random sample, Gamma prior distribution $\longrightarrow$ Gamma posterior distribution \item Normal random sample (known variance), Normal prior distribution $\longrightarrow$ Normal posterior distribution \end{itemize} In these examples, the prior distribution and the posterior distribution come from the same family. This leads us to the following definition. \begin{definition}~\\ Suppose that data $\xvec$ are to be observed with distribution $f(\xvec|\theta)$. A family $\mathfrak{F}$ of prior distributions for $\theta$ is said to be {\it conjugate} to $f(\xvec|\theta)$ if for every prior distribution $\pi(\theta)\in\mathfrak{F}$, the posterior distribution $\pi(\theta|\xvec)$ is also in $\mathfrak{F}$. \end{definition} Notice that the conjugate family depends crucially on the model chosen for the data $\xvec$. For example, the only family conjugate to the model ``random sample from a Poisson distribution'' is the Gamma family. \subsection{Vague prior knowledge/prior ignorance} If we have very little or no prior information about the model parameters $\theta$, we must still choose a prior distribution in order to operate Bayes Theorem. Obviously, it would be sensible to choose a prior distribution which is not concentrated about any particular value, that is, one with a very large variance. In particular, most of the information about $\theta$ will be passed through to the posterior distribution via the data, and so we have $\pi(\theta|\xvec)\sim f(\xvec|\theta)$. We often represent vague prior knowledge by using a prior distribution which is conjugate to the model for $\xvec$ and which is as diffuse as possible, that is, has as large a variance as possible. We may attempt to represent prior ignorance using a \textit{Jeffreys prior}, specified in terms of Fisher's Information: \begin{equation} \label{eq:fisher} I(\theta)=E_{\Xvec|\theta}\left[- \frac{\partial^2}{\partial\theta^{2}}\, \log f(\Xvec|\theta)\right]. \end{equation} Specifically, Jeffreys recommended that we represent prior ignorance by the prior distribution \begin{equation} \pi(\theta)\propto \sqrt{I(\theta)}. \end{equation} Such a prior distribution is known as a Jeffreys prior distribution. \begin{example}~\\ In constructing the prior in Example 1.1, did we assume (i) substantial prior knowledge; (ii) limited prior knowledge; or (iii) prior ignorance? \framebox[17cm]{ \parbox[t]{16cm}{ \vspace{0.4cm} \huge \ding{46} \normalsize \textcolor{blue}{ Limited prior knowledge -- for mathematical convenience we used a conjugate prior, but we obtained the parameters in this prior through prior information collected from the epidemiologist. \vspace{1.5cm} } }} \end{example} \begin{example}~\\ Suppose the parameter $\theta$ represents a probability. What might be the simplest prior to use here, if we are prior ignorant? \framebox[17cm]{ \parbox[t]{16cm}{ \vspace{0.4cm} \huge \ding{46} \normalsize \textcolor{blue}{ Rather than compute a Jeffreys prior, we might simply choose a Uniform $U(0, 1)$ prior here; that is, $\theta \sim U(0,1)$, with \[ \pi(\theta) = 1. \] This would give all outcomes across the permissible range of values for $\theta$ equal density.} \vspace{1cm} }} \end{example} \newpage \begin{example}~\\ Suppose we have a random sample from a $N(\mu,1/\tau)$ distribution (with $\tau$ known). Determine the posterior distribution assuming a vague prior for $\mu$. \framebox[17cm]{ \parbox[t]{16cm}{ \vspace{0.4cm} \huge \ding{46} \normalsize \textcolor{blue}{ The conjugate prior distribution is a normal distribution. We have already seen that if the prior is $\mu\sim N(b,1/d)$ then the posterior distribution is $\mu|\xvec\sim N(B,1/D)$ where \begin{equation*} B=\frac{db+n\tau\bar x}{d+n\tau}\qquad\text{and}\qquad D=d+n\tau. \end{equation*} If we now make our prior knowledge vague about $\mu$ by letting the prior variance tend to infinity ($d\to 0$), we obtain \[ B\to\bar x\qquad\text{and}\qquad D\to n\tau. \] Therefore, assuming vague prior knowledge for $\mu$ results in a $N\{\bar x,1/(n\tau)\}$ posterior distribution. Notice that the posterior mean is the sample mean (the likelihood mode) and that the posterior variance $1/(n\tau)\to 0$ as $n\to\infty$.} \vspace{1.5cm} }} \end{example} \begin{example}~\\ Suppose we have a random sample from a Poisson distribution, that \label{ex:vaguegamma} is, $X_i|\theta\sim Po(\theta)$, $i=1,2,\ldots,n$ (independent). Determine the posterior distribution assuming a vague prior for $\theta$. \framebox[17cm]{ \parbox[t]{16cm}{ \vspace{0.4cm} \huge \ding{46} \normalsize \textcolor{blue}{ The conjugate prior distribution is a Gamma distribution. Recall that a $Ga(g,h)$ distribution has mean $m=g/h$ and variance $v=g/h^2$. Rearranging these formulae we obtain \begin{equation*} g=\frac{m^2}{v}\quad\quad\text{and}\quad\quad h=\frac{m}{v}. \end{equation*} Clearly $g\to 0$ and $h\to 0$ as $v\to\infty$ (for fixed~$m$). We have seen how taking a $Ga(g,h)$ prior distribution results in a $Ga(g+n\bar x,h+n)$ posterior distribution. Therefore, taking a vague prior distribution will give a $Ga(n\bar x,n)$ posterior distribution. Note that the posterior mean is $\bar x$ (the likelihood mode) and that the posterior variance $\bar x/n\to 0$ and $n\to\infty$.} \vspace{2cm}}} \end{example} \begin{example}~\\ Suppose we have a random sample from an exponential distribution, that \label{ex:38} is, $X_i|\theta\sim \textrm{\it Exp}(\theta)$, $i=1,2,\ldots,n$ (independent). Determine the Jeffreys prior for this model. \framebox[17cm]{ \parbox[t]{16cm}{ \vspace{0.4cm} \huge \ding{46}\normalsize \textcolor{blue}{ Recall that $%\[ L(\theta|\xvec)=f_{\Xvec}(\xvec|\theta)=\theta^n e^{-n\bar x\theta}, $ %\] and therefore \begin{align*} \log &L(\theta|\xvec)=n\log\theta -n\bar x\theta \\ &\Rightarrow\quad \frac{\partial}{\partial\theta} \log L(\theta|\xvec)=\frac{n}{\theta}-n\bar x \\ &\Rightarrow\quad \frac{\partial^2}{\partial\theta^2} \log L(\theta|\xvec)=-\frac{n}{\theta^2} \\ &\Rightarrow\quad I(\theta)=E_{\Xvec|\theta} \left[-\frac{\partial^2}{\partial\theta^2} \log L(\theta|\Xvec)\right]=\frac{n}{\theta^2}. \end{align*} Hence, the Jeffreys prior for this model is \begin{align*} \pi(\theta)&\propto I(\theta)^{1/2} \\ &\propto \frac{\sqrt{n}}{\theta},\quad\quad\theta>0 \\ &\propto \frac{1}{\theta},\quad\quad\theta>0. \end{align*}} \vspace{6cm} }} Notice that this distribution is improper since $\int_0^\infty d\theta/\theta$ is a divergent integral, and so we cannot find a constant which ensures that the density function integrates to one. Notice also that this density is, in fact, a limiting form of a $\textrm{\it Ga}(g,h)$ density (ignoring the integration constant) since \begin{equation*} \frac{h^g\,\theta^{g-1}e^{-h\theta}}{\Gamma(g)}\propto \theta^{g-1}e^{-h\theta} \to \frac{1}{\theta},\quad\quad\text{as }g\to 0,~h\to 0. \end{equation*} Therefore, we obtain the same posterior distribution whether we adopt the Jeffreys prior or vague prior knowledge. \end{example} \section{The asymptotic posterior distribution} If we have a statistical model $f(\xvec|\theta)$ for data $\xvece$, together with a prior distribution $\pi(\theta)$ for $\theta$ then \begin{align*} &\sqrt{J(\thetahat)}~(\theta-\thetahat)|\xvec \stackrel{\cal D} \longrightarrow N(0,1) \quad\quad\text{as }n\to\infty, \\ \intertext{where $\thetahat$ is the likelihood mode and $J(\theta)$ is the observed information} &J(\theta)=-\frac{\partial^2}{\partial\theta^2}\, \log f(\xvec|\theta). \end{align*} This means that, with increasing amounts of data, the posterior distribution looks more and more like a normal distribution. The result also gives us a useful approximation to the posterior distribution for $\theta$ when $n$ is large: \[ \theta|\xvec\sim N\{\hat\theta,J(\hat\theta)^{-1}\} \quad\quad\text{approximately}. \] Note that this limiting result is similar to one used in Frequentist statistics for the distribution of the maximum likelihood estimator, namely \begin{align*} &\sqrt{I(\theta)}~(\thetahat-\theta)\stackrel{\cal D} \longrightarrow N(0,1) \quad\quad\text{as }n\to\infty, \end{align*} where Fisher's information $I(\theta)$ is the expected value of the observed information, where the expectation is taken over the distribution of $\Xvec|\theta$, that is, $I(\theta)=E_{\Xvec|\theta}[J(\theta)]$. You may also have seen this result written as an approximation to the distribution of the maximum likelihood estimator in large samples, namely \[ \thetahat\sim N\{\theta,I(\theta)^{-1}\} \quad\quad\text{approximately}. \] \begin{example}~\\ Suppose we have a random sample from a $N(\mu,1/\tau)$ distribution (with $\tau$ known). Determine the asymptotic posterior distribution for~$\mu$. Recall that \[ f(\xvec|\mu)= \left(\frac{\tau}{2\pi}\right)^{n/2} \exp\left\{-\frac{\tau}{2}\sum_{i=1}^n (x_i-\mu)^2\right\}, \] and therefore \begin{align*} \log\,&f(\xvec|\mu)=\frac{n}{2}\log\tau-\frac{n}{2}\log(2\pi) -\frac{\tau}{2}\sum_{i=1}^n (x_i-\mu)^2 \\ &\Rightarrow\quad \frac{\partial}{\partial\mu} \log f(\xvec|\mu)=-\frac{\tau}{2}\times\sum_{i=1}^n -2(x_i-\mu) =\tau\sum_{i=1}^n (x_i-\mu) =n\tau(\bar x-\mu)\\ &\Rightarrow\quad \frac{\partial^2}{\partial\mu^2} \log f(\xvec|\mu)=-n\tau \qquad\Rightarrow\quad J(\mu)=-\frac{\partial^2}{\partial\mu^2} \log f(\xvec|\mu)=n\tau. \end{align*} \framebox[17cm]{ \parbox[t]{16cm}{ \vspace{0.4cm} \huge \ding{46}\normalsize \textit{... Solution to Example 1.10} \textcolor{blue}{ We have \begin{align*} \frac{\partial}{\partial\mu} \log f(\xvec|\mu)=0 \quad\quad&\Longrightarrow\quad\quad \hat\mu=\bar x \\ &\Longrightarrow\quad\quad J(\hat\mu)=n\tau \\ &\Longrightarrow\quad\quad J(\hat\mu)^{-1}=\frac{1}{n\tau}. \end{align*} Therefore, for large $n$, the (approximate) posterior distribution for $\mu$ is \[ \mu|\xvec\sim N\left(\bar x,\,\frac{1}{n\tau}\right). \] Here the asymptotic posterior distribution is the same as the posterior distribution under vague prior knowledge.} \vspace{2.3cm} }} \end{example} \section{Bayesian inference} The posterior distribution $\pi(\theta|\xvec)$ summarises all our information about $\theta$ to date. However, sometimes it is helpful to reduce this distribution to a few key summary measures. \subsection{Estimation} \label{sec:univariate-estimation} \subsubsection*{Point estimates} There are many useful summaries for a typical value of a random variable with a particular distribution; for example, the mean, mode and median. The mode is used more often as a summary than is the case in frequentist statistics. \subsubsection*{Confidence intervals/regions} A more useful summary of the posterior distribution is one which also reflects its variation. For example, a $100(1-\alpha)\%$ {\it Bayesian confidence interval} for $\theta$ is any region $C_\alpha$ that satisfies $Pr(\theta\in C_\alpha|\xvec)=1-\alpha$. If $\theta$ is a continuous quantity with posterior probability density function $\pi(\theta|\xvec)$ then \begin{equation*} \int_{C_\alpha} \pi(\theta|\xvec)\,d\theta = 1-\alpha. \end{equation*} The usual correction is made for discrete $\theta$, that is, we take the largest region $C_\alpha$ such that $Pr(\theta\in C_\alpha|\xvec)\leq 1-\alpha$. Bayesian confidence intervals are sometimes called {\it credible regions} or {\it plausible regions}. Clearly these intervals are not unique, since there will be many intervals with the correct probability coverage for a given posterior distribution. A $100(1-\alpha)\%$ \textit{highest density interval} (HDI) for $\theta$ is the region \[ C_\alpha=\{\theta:~\pi(\theta|\xvec)\geq\gamma\} \] where $\gamma$ is chosen so that $Pr(\theta\in C_\alpha|\xvec)=1-\alpha$. This region is sometimes called a \textit{most plausible Bayesian confidence interval}. If the posterior distribution has many modes then it is possible that the HDI will be the union of several disjoint regions. Also, if the posterior distribution is unimodal (has one mode) and symmetric about its mean then the HDI is an equi-tailed interval, that is, takes the form $C_\alpha=(a,b)$, where $Pr(\thetab|\xvec)=\alpha/2$; see Figure~\ref{fig:hdi}. \begin{figure}[!h] \begin{center} \includegraphics[width=0.8\textwidth]{hdi} \caption{Construction of an HDI for a symmetric posterior density} \label{fig:hdi} \end{center} \end{figure} %%Hack for vspace %\newpage \subsubsection*{Interpretation of confidence intervals/regions} Suppose $C_B$ is a 95\% Bayesian confidence interval for $\theta$ and $C_F$ is a 95\% frequentist confidence interval for $\theta$. These intervals do not have the same interpretation: \begin{myitemize} \item the probability that $C_B$ contains $\theta$ is 0.95; \item the probability that $C_F$ contains $\theta$ is either 0 or 1 --- since $\theta$ does not have a (non-degenerate) probability distribution; \item the interval $C_F$ covers the true value $\theta$ on 95\% of occasions --- in repeated applications of the formula. \end{myitemize} \begin{example}~\\ Suppose we have a random sample $\xvece$ from a $N(\mu,1/\tau)$ distribution (where $\tau$ is known). We have seen that, assuming vague prior knowledge, the posterior distribution is $\mu|\xvec\sim N\{\bar x,1/(n\tau)\}$. Determine the $100(1-\alpha)\%$ HDI for~$\mu$. \framebox[17cm]{ \parbox[t]{16cm}{ \vspace{0.4cm} \huge \ding{46}\normalsize \textcolor{blue}{ This distribution has a symmetric bell shape and so the HDI is an equi-tailed interval $C_\alpha=(a,b)$ with $Pr(\mub|\xvec)=\alpha/2$, that is, \[ a=\bar x - \frac{z_{\alpha/2}}{\sqrt{n\tau}} \quad\quad\text{and}\quad\quad b=\bar x + \frac{z_{\alpha/2}}{\sqrt{n\tau}}, \] where $z_\alpha$ is the upper $\alpha$-quantile of the $N(0,1)$ distribution. For example, the 95\% HDI for~$\mu$ is \[ \left(\bar x - \frac{1.96}{\sqrt{n\tau}},\, \bar x + \frac{1.96}{\sqrt{n\tau}}\right). \] } }} Note that this interval is numerically identical to the 95\% frequentist confidence interval for the (population) mean of a normal random sample with known variance. However, the interpretation is very different. \end{example} \begin{example}~\\ Recall Example~\ref{ex:botulism} on the number of cases of foodbourne botulism in England and Wales. The data were modelled as a random sample from a Poisson distribution with mean~$\theta$. Using a $Ga(2,1)$ prior distribution, we found the posterior distribution to be $\theta|\xvec\sim Ga(8,9)$ (see Figure~\ref{fig:botulismpost2}). Determine the $100(1-\alpha)\%$ HDI for~$\theta$. \begin{figure}[ht] \begin{center} \includegraphics[width=0.6675\textwidth]{botulism2} \caption{Posterior density for $\theta$} \label{fig:botulismpost2} \end{center} \end{figure} \newpage \framebox[17cm]{ \parbox[t]{16cm}{ \vspace{0.4cm} \huge \ding{46}\normalsize \textit{... Solution to Example 1.12} \textcolor{blue}{ The HDI must take the form $C_\alpha=(a,b)$ if it is to include the values of $\theta$ with the highest probability density. Suppose that $F(\cdot)$ and $f(\cdot)$ are the posterior distribution and density functions. Then the end-points $a$ and $b$ must satisfy \begin{align*} Pr(a<\theta