\documentclass[12pt,a4paper]{article} \usepackage{amsthm,amsmath,amsfonts} \usepackage{latexsym} \usepackage{fancyhdr} \usepackage{enumerate} %\usepackage{marginnote} \usepackage{graphicx} \usepackage{skull} \usepackage{hyperref} \usepackage{fancybox} \usepackage{xspace} %%From web for code highlighting:: \usepackage{listings} \usepackage{cmbright} \addtolength{\parskip}{1ex} \baselineskip = 6mm %\parskip = 2mm \textheight 23.6cm \textwidth 16.6cm \oddsidemargin 0.0cm \evensidemargin 0.0cm \topmargin -1.4cm %\headsep -.60in \setlength{\parindent}{0cm} \pagestyle{fancy} \lhead{} \chead{} \rhead{} \lfoot{} \cfoot{\thepage} \rfoot{} \rhead{} \newcommand{\pr}[2][]{\ensuremath{\mathrm{Pr}_{#1}\!\left(#2\right)}} \newenvironment{myhighlight}{\textit}{} \renewenvironment{emph}{\textbf}{} \newcommand{\Exp}[2][]{\ensuremath{\mathrm{Exp}\left(#2\right)}} \newcommand{\Bin}[2][]{\ensuremath{\mathrm{Bin}\left(#2\right)}} \newcommand{\Po}[2][]{\ensuremath{\mathrm{Po}\left(#2\right)}} \newcommand{\cov}[2][]{\ensuremath{\mathrm{Cov}\left(#2\right)}} \newcommand{\cor}[2][]{\ensuremath{\mathrm{Cor}\left(#2\right)}} \newcommand{\diag}[2][]{\ensuremath{\mathrm{diag}\left(#2\right)}} \newcommand{\given}{\ensuremath{\ |\ }} \newcommand{\stationary}{\ensuremath{\vec{\pi}}} \newcommand{\identity}{\ensuremath{\vec{I}}} \renewcommand{\vec}[1]{\underline{#1}} \newcommand{\remove}[1]{} %\newcommand{\cover}[2]{\hi{#1}} \newcommand{\cover}[2]{#2} \newcommand{\hi}[2][]{\textbf{#2}} \newcommand{\bigfont}{\fontsize{1.45cm}{1.5cm}\selectfont} \newcommand{\announce}{\setbeamercolor{structure}{fg=green!60!black}} \newcommand{\sol}{\setbeamercolor{structure}{fg=purple!80!black}} \newcommand{\normal}{\setbeamercolor{structure}{fg=blue!60!green}} \newcommand{\vsp}{\vspace{0.3cm}} \newcommand{\vspl}{\vspace{0.15cm}} \newcommand{\xvec}{\boldsymbol{x}} \newcommand{\yvec}{\vec{y}} \newcommand{\Xvec}{\boldsymbol{X}} \newcommand{\Yvec}{\vec{Y}} \newcommand{\Evec}{\vec{E}} \newcommand{\tvec}{\bodsymbol{t}} \newcommand{\Tvec}{\bodsymbol{T}} \newcommand{\mutau}{\begin{pmatrix}\mu\\\tau\end{pmatrix}} \newcommand{\thetavec}{\vec{\theta}} \newcommand{\thetavechat}{\vec{\hat\theta}} \newcommand{\thetahat}{\hat\theta} \newcommand{\xvece}{\boldsymbol{x}=(x_1,x_2,\ldots,x_n)^T} \newcommand{\Xvece}{\boldsymbol{X}=(X_1,X_2,\ldots,X_n)^T} \newcommand{\thetavece}{\vec{\theta}=(\theta_1,\theta_2,\ldots,\theta_p)^T} \newcommand{\exxvec}{$x_1,x_2,\ldots,x_n$\xspace} \newcommand{\exxvecdot}{$x_1,x_2,\ldots,x_n$.\ } \newcommand{\exXvec}{$X_1,X_2,\ldots,X_n$\ } \newcommand{\muvec}{\vec{\mu}} \newtheoremstyle{exmp}% hnamei {3pt}% hSpace abovei {3pt}% hSpace belowi {}% hBody fonti {}% hIndent amounti1 {\bf}% hTheorem head fonti {:}% hPunctuation after theorem headi {.5em}% hSpace after theorem headi2 {}% hTheorem head spec (can be left empty, meaning `normal')i \theoremstyle{exmp} \newtheorem{example}{\textcolor{purple}{\textbf{Example}}}[section] %\def \var{{\textrm{Var}}} \newcommand{\E}[2][]{\ensuremath{\mathrm{E}\left[#2\right]}} \def \R{\mathbb{R}} \def \Z{\mathbb{Z}} \def \proof{\noindent{\bf Proof: }} \def \eproof{\qed\vskip .3cm} \def \union{\cup} \def \intersect{\cap} \newcommand{\ex}[2][]{\ensuremath{\mathrm{E}_{#1}\!\left[#2\right]}} \newcommand{\var}[2][]{\ensuremath{\mathrm{Var}_{#1}\left(#2\right)}} %\renewcommand{\det}[1][]{\ensuremath{\mathrm{det}\left(#1\right)}} \newtheorem{defn}{Definition} %\newtheorem*{example}{Example} \newtheorem{result}{Result}[section] \newtheorem{algorithm}{Algorithm}[section] \usepackage{Sweave} \usepackage[usenames, dvipsnames]{xcolor} \begin{document} %%From web for code/highlighting \definecolor{mygreen}{rgb}{0,0.6,0} \definecolor{mygray}{rgb}{0.5,0.5,0.5} \definecolor{mymauve}{rgb}{0.58,0,0.82} \definecolor{light-gray}{gray}{0.95} \definecolor{LEEgreen}{rgb}{0.04,0.667,0.137} \lstset{ backgroundcolor=\color{light-gray}, % choose the background color; you must add \usepackage{color} or \usepackage{xcolor}; should come as last argument basicstyle=\footnotesize, % the size of the fonts that are used for the code breakatwhitespace=false, % sets if automatic breaks should only happen at whitespace breaklines=true, % sets automatic line breaking captionpos=b, % sets the caption-position to bottom commentstyle=\color{mygreen}, % comment style deletekeywords={...}, % if you want to delete keywords from the given language escapeinside={\%*}{*)}, % if you want to add LaTeX within your code extendedchars=false, % lets you use non-ASCII characters; for 8-bits encodings only, does not work with UTF-8 frame=single, % adds a frame around the code keepspaces=true, % keeps spaces in text, useful for keeping indentation of code (possibly needs columns=flexible) keywordstyle=\color{blue}, % keyword style language=R, % the language of the code morekeywords={*,...}, % if you want to add more keywords to the set numbers=left, % where to put the line-numbers; possible values are (none, left, right) numbersep=5pt, % how far the line-numbers are from the code numberstyle=\tiny\color{mygray}, % the style that is used for the line-numbers rulecolor=\color{black}, % if not set, the frame-color may be changed on line-breaks within not-black text (e.g. comments (green here)) showspaces=false, % show spaces everywhere adding particular underscores; it overrides 'showstringspaces' showstringspaces=false, % underline spaces within strings only showtabs=false, % show tabs within strings adding particular underscores stepnumber=1, % the step between two line-numbers. If it's 1, each line will be numbered stringstyle=\color{mymauve}, % string literal style tabsize=2, % sets default tabsize to 2 spaces title=\lstname % show the filename of files included with \lstinputlisting; also try caption instead of title } %\SweaveOpts{concordance=TRUE} \pagenumbering{gobble} \setcounter{page}{-3} % to get contents page labelled (i) %%%Practicals %\section*{Aims} %This introductory practical is split into two parts: %\begin{itemize} %\item In the first hour or so, we will focus on the basics of R, revising some of the key concepts from last year's MAS1802 material %\item In the second hour we will focus on material from this week's lecture -- using R to simulate from standard probability distributions, as well as more complex models including mixture distributions %\end{itemize} %\vspace{0.25cm} %\noindent Please read through all of this handout \textit{very carefully}, and follow the instructions \textit{meticulously}. Whenever %you are typing commands into R, make sure you do so \textit{extremely carefully}. If you get an error message, %try typing the code in again, just to make sure you've not made a spelling mistake or other simple typo. A %common error is to forget that R is case sensitive, so use upper/lower cases exactly when you're told to %do so. %\vspace{0.25cm} %\noindent For all R practicals, make sure you've got your lecture notes handy. %\section*{Getting started} %\noindent Log-in to a PC in the usual way. Now open Rstudio. If you can't remember how to do this, click on the ``Start'' button in the bottom %left-hand-corner and search for ``rstudio''. Once you've opened Rstudio you should see four panels; if there is no text editor (usually the top-left panel), click on the green button with a white plus symbol and open a new R script. Remember -- to execute a line of code in the text editor, hit \texttt{Ctrl}-\texttt{Enter}. %\vspace{0.25cm} %\noindent It is good practice to always use the text editor, rather than type into the console window (bottom- %left panel) directly. This is especially true in future practicals when you may be building large functions; it %is much easier to de-bug and edit code in the text editor than in the console window. %\vspace{0.25cm} %\noindent \textbf{Hint: To commit several lines of code in one go, highlight your block of code in the text editor and hit \texttt{Ctrl}-\texttt{Enter} just once!} %\vspace{0.25cm} %\noindent Now you've got Rstudio open, you're ready to go! Try out \textit{all} of the following exercises. %\newpage \lhead{MAS3902} \rhead{Bayesian Inference} \section*{Bayesian Inference: Assignment 2} \noindent Answers to the starred (\textbf{*}) questions should be submitted to Canvas by no later than \textcolor{red}{\bf 4pm, Friday 8th March 2024}. Solutions to non-assessed questions will be made available in class and/or in Canvas. \subsection*{Question 1*} Re-visit question 1 from assignment 1. Recall that we were interested in $\theta$, the probability that particulate matter in the air as a result of a dust storm, results in an environmental health warning being issued. \begin{itemize} \item [(a)] This question concerns Bayesian confidence intervals: \begin{itemize} \item [(i)] Assuming the informative prior used in assignment 1 (not the ignorance prior), which of the following represent valid 95\% posterior Bayesian confidence intervals for $\theta$? Explain your answer. \textcolor{gray}{\bf [3 marks]} \[ (0.2463, 0.5331) \qquad (0.1342, 0.4359) \qquad (0.1102, 0.5021) \qquad (0.1422, 0.4461) \] \item [(ii)] Which interval listed in part (a) is a 95\% posterior highest density interval for $\theta$? Explain your answer. \textcolor{gray}{\bf [2 marks]} \item [(iii)] Find the 95\% prior highest density interval for $\theta$, and compare it to that obtained from the posterior in part (ii). \textcolor{gray}{\bf [2 marks]} \end{itemize} \item [(b)] This question is about the Bayesian posterior predictive distribution: \begin{itemize} \item [(i)] For the next ten major dust storms to hit Greece and Cyprus this year, find the predictive probability distribution for the number $Y$ that will result in an environmental health warning. \textcolor{gray}{\bf [9 marks]} \textcolor{gray}{ Hint: The random variable $X$ follows a Beta-binomial $BetaBin(n, a, b)$ distribution ($n$ positive integer, $a>0, b>0$) if it has probability function \[ f(x|n, a, b) = \left(\begin{array}{c} n\\ x \end{array} \right) \frac{B(x+a, b+n-x)}{B(a,b)}, \quad x=0, 1, \ldots, n, \] where $B(a, b)$ is the beta function and \[ B(a, b) = \frac{(a-1)!(b-1)!}{(a+b-1)!}, \quad \text{for $a$ and $b$ positive integers.} \]} \item [(ii)] Use your answer to part (i) to find a 95\% prediction set for the number of dust storms (of the next ten observed) to trigger an environmental health warning. \textcolor{gray}{\bf [4 marks]} \end{itemize} \end{itemize} \subsection*{Question 2} The negative binomial distribution is used to model scenarios in which we observe the number of independent (success-fail) trials needed before we see a successful trial; the trials must have the same success probability. Suppose \exxvec are a random sample from a negative binomial $NegBin(k,\theta)$ distribution, where $k$ is known. \begin{itemize} \item[(i)] Verify that the congugate prior distribution is a Beta distribution. \item[(ii)] Determine the choice of parameters $g$ and $h$ for the $Beta(g,h)$ distribution that give it maximal variance.\\ Hint: reparameterise the distribution in terms of its mean $m=g/(g+h)$ and $s=g+h$; determine the choice of $m$ and $s$ that maximises the variance of the distribution, and hence the choice of $g$ and $h$. \end{itemize} Determine the posterior distribution for $\theta$ assuming \begin{itemize} \item[(iii)] vague prior knowledge; \item[(iv)] a very large sample. \end{itemize} \subsection*{Question 3} The Rayleigh distribution is often used to measure variability in magnetic resonance imaging (MRI). Suppose \exxvec are a random sample from a Rayleigh $R(\theta)$ distribution. Determine the posterior distribution for $\theta$ assuming \begin{enumerate}[(i)] \item vague prior knowledge; \item a very large sample. \end{enumerate} \subsection*{Question 4} There were 4 ``significant'' flooding events in the UK in 2022. The posterior predictive distribution for $Y$, the number of significant flooding events in 2023, is shown below: \vspace{-0.5cm} \begin{small} \begin{center} \begin{tabular}{|c|c|c|c|c|c|c|c|c|c|c|c|}\hline $y$ &0 &1 &2 &3 &4 &5 &6 &7 &8 &9 &10 \\ \hline $f(y|x=4)$ &0.00 &0.01&0.04&0.10 &0.13 &0.45 &0.17 &0.05 &0.03 &0.01 &0.01 \\ \hline \end{tabular} \end{center} \end{small} \begin{itemize} \item [(a)] Find the most plausible 90\% prediction set for $Y$. \item [(b)] Find $E(y|x=4)$ and $SD(y|x=4)$. \end{itemize} \newpage \section*{Solutions} \subsection*{Question 1} \begin{itemize} \item [(a)] Bayesian confidence intervals: \begin{itemize} \item [(i)] Recall from question 1 in assignment 1 that the posterior for $\theta$ is $\theta|\boldsymbol{x} \sim Be(9,23)$. A valid 95\% Bayesian confidence interval $(a, b)$ for $\theta|\boldsymbol{x}$ would be such that $Pr(a<\theta pbeta(0.5331, 9, 23) - dbeta(0.2463, 9, 23) [1] 0.6501951 > pbeta(0.4359, 9, 23) - dbeta(0.1342, 9, 23) [1] 0.9499906 > pbeta(0.5021, 9, 23) - dbeta(0.1102, 9, 23) [1] 0.9899973 > pbeta(0.4461, 9, 23) - dbeta(0.1422, 9, 23) [1] 0.9500078 \end{verbatim} Thus, both $(0.1342, 0.4359)$ and $(0.1422, 0.4461)$ are valid 95\% Bayesian confidence intervals. \textcolor{red}{[1 mark for correct identification, 2 marks for reasoning]} \item [(ii)] We know from the course that a highest density interval also has equal density ($\gamma$) at both $a$ and $b$. Again, we can check in R (or we could do so by hand!). For the first interval: \begin{verbatim} > dbeta(0.1342, 9, 23) [1] 0.8015702 > dbeta(0.4359, 9, 23) [1] 0.8010577 \end{verbatim} And for the second interval: \begin{verbatim} > dbeta(0.1422, 9, 23) [1] 1.038543 > dbeta(0.4461, 9, 23) [1] 0.6451823 \end{verbatim} Thus, only the first interval $(0.1342, 0.4359)$ satisfies $\pi(\theta=a|\boldsymbol{x})=\pi(\theta=b|\boldsymbol{x})=\gamma$, and so $(0.1342, 0.4359)$ is the highest density interval. \textcolor{red}{[2 marks]} We could check this using the \texttt{nclbayes} package: \begin{verbatim} > hdiBeta(0.95, 9, 23) [1] 0.1341826 0.4358943 \end{verbatim} \item [(iii)] We could use the \texttt{nclbayes} package: \begin{verbatim} > hdiBeta(0.95, 5, 17) [1] 0.06921858 0.39949020 \end{verbatim} Thus, the prior HDI is $(0.069, 0.399)$ \textcolor{red}{[1 mark]}-- similar to that for the posterior (it overlaps!), but shifted to the left -- suggesting we are now slightly more inclined to believe that a dust storm will trigger a health warning. \textcolor{red}{[1 mark]} \end{itemize} %\newpage \item [(b)] Predictive distribution: \begin{itemize} \item [(i)] Suppose the number of dust storms triggering a warning in the future is $Y|\theta\sim Bin(10,\theta)$. The predictive probability function of $Y$ is, for $y=0,1,2,\ldots,10,$ \begin{eqnarray*} f(y|x=4) &=&\int_\Theta f(y|\theta)\,\pi(\theta|x=4)\,d\theta \\ & & \\ &=&\int_0^1 \binom{10}{y}\theta^y(1-\theta)^{10-y} \times \frac{\theta^{8}(1-\theta)^{22}}{B(9,23)}\,d\theta \\ & & \\ &=&\binom{10}{y}\frac{1}{B(9,23)}\int_0^1 \theta^{y+8}(1-\theta)^{32-y}\,d\theta \quad \text{\textcolor{red}{[2 marks]}}\\ & & \\ &=& \binom{10}{y}\frac{\color{blue}B(y+9, 33-y)}{B(9,23)} \int_0^1 \frac{\theta^{(y+9)-1}(1-\theta)^{(33-y)-1}}{{\color{blue}B(y+9, 33-y)}}\,d\theta \quad \text{\textcolor{red}{[4 marks]} }\\ &=&\binom{10}{y}\frac{B(y+9,33-y)}{B(9,23)}. \quad \text{\textcolor{red}{[1 mark]}} \end{eqnarray*} That is, $Y|x=4\sim BetaBin(10,9,23)$. \textcolor{red}{[2 mark]} \item [(ii)] Probably best to write an R function? \begin{verbatim} > betabin = function(y, n, a, b){ numerator = beta(y + a, n + b -y) denominator = beta(a, b) answer = choose(10,y)*numerator/denominator return(answer)} > cbind(0:10,round(betabin(0:10, 10, 9, 23),3)) [,1] [,2] [1,] 0 0.058 [2,] 1 0.162 [3,] 2 0.235 [4,] 3 0.230 [5,] 4 0.166 [6,] 5 0.093 [7,] 6 0.040 [8,] 7 0.013 [9,] 8 0.003 [10,] 9 0.000 [11,] 10 0.000 \end{verbatim} The above is our predictive distribution, that is, $f(y|x=4)$. For the prediction set, we need to include values for $y$ with greatest predictive probability, such that the sum of their probabilities is close to 0.95. We include, in order (and with cumulative probabilities in parentheses): $[(2(0.235), 3(0.465), 4(0.631), 1(0.793), 5(0.886), 0(0.944)]$, i.e. $[0, 5]$. \textcolor{red}{[3 marks for method, 1 mark for final answer of $[0, 5]$.} \end{itemize} \end{itemize} %\newpage \subsection*{Question 2} \begin{itemize} \item [(i) ]If the prior distribution is $\theta\sim Beta(g,h)$ then, using Bayes Theorem, the posterior density is \begin{align*} \pi(\theta|\xvec)&\propto\pi(\theta)\,f(\xvec|\theta) \\ &\propto \theta^{g-1}(1-\theta)^{h-1}\times \prod_{i=1}^n\theta^k(1-\theta)^{x_i-k}, \quad\quad0<\theta<1 \\ &\propto \theta^{g+nk-1}(1-\theta)^{h+n\bar x-nk-1}, \quad\quad0<\theta<1 \\ \noalign{\medskip} i.e.&~\theta|\xvec\sim Beta(G=g+nk,H=h+n\bar x-nk). \end{align*} Therefore, both the prior and posterior distributions are Beta distributions and so the Beta is the conjugate prior distribution. \item[(ii)] As the beta distribution restricts values to the range (0,1), there is a finite upper limit to the variance. If $\theta\sim Beta(g,h)$ then it has variance $Var(\theta)=gh/\{(g+h)^2(g+h+1)\}$. Note that this distribution requires its parameters $g>0$ and $h>0$. If we now transform the parameters to $m=g/(g+h)$ and $s=g+h$ then $Var(\theta)=m(1-m)/(s+1)$. Note that, as this is a 1--1 transformation, we have a free choice of values for $m$ and $s$ subject to the restrictions $00$. Clearly, this ratio is maximised when the numerator is maximised and the denominator is minimised. Now $m(1-m)$ is maximised at $m=1/2$, and $(s+1)$ is minimised by taking $s$ as small as possible. In terms of $g$ and $h$, these conditions are $g=h$ and $g+h\to 0$, that is, $g=h\to 0$. Intuitively, the maximum variance is achieved when the probability density is pushed to the extremes of the range, that is, equal mass at $\theta=0$ and $\theta=1$ -- this distribution is obtained in the limit $g=h\to 0$. \item[(iii)] We make the Beta prior distribution vague by making it have as large a variance as possible. This is achieved, from (ii) above, by letting $g=h\to 0$. Hence, using the result in (i), the posterior distribution assuming vague prior information is \begin{equation*} \theta|\xvec\sim Beta(nk,n\bar x-nk). \end{equation*} \item[(iv)] The asymptotic posterior distribution (as $n\to\infty$) is \begin{equation*} \theta|\xvec\sim N\left(\hat\theta,J(\hat\theta)^{-1}\right), \end{equation*} where \begin{equation*} J(\theta)=-\frac{\partial^2}{\partial\theta^2} \log f(\xvec|\theta)=\frac{nk}{\theta^2}+\frac{n(\bar x-k)}{(1-\theta)^2}. \end{equation*} Now \begin{eqnarray*} \frac{\partial}{\partial\theta} \log f(\xvec|\theta)=0 \quad\quad&\Longrightarrow&\quad\quad \frac{nk}{\hat\theta}-\frac{n(\bar x-k)}{1-\hat\theta}=0 \\ &\Longrightarrow&\quad\quad \hat\theta=\frac{k}{\bar x} \\ &\Longrightarrow&\quad\quad J(\hat\theta)=\frac{n\bar x^3}{k(\bar x-k)} \\ &\Longrightarrow&\quad\quad J(\hat\theta)^{-1}=\frac{k(\bar x-k)}{n\bar x^3}. \end{eqnarray*} Therefore, for large $n$, the posterior distribution for $\theta$ is \begin{equation*} \theta|\xvec\sim N\left(\frac{k}{\bar x},\,\frac{k(\bar x-k)}{n\bar x^3}\right) \quad\quad\text{approximately}. \end{equation*} \end{itemize} \subsection*{Question 3} \begin{itemize} \item [(i)] The likelihood function is \[ f(\xvec|\theta)=\prod_{i=1}^n 2x_i\theta e^{-\theta x_i^2} \propto\theta^n\exp\left\{-\theta\sum_{i=1}^n x_i^2\right\} =\theta^n\exp\left\{-n\overline{x^2}\theta\right\}, \] where $\overline{x^2}=\sum_{i=1}^n x_i^2/n$ is the mean of the squared observations. Therefore the conjugate prior distribution is a Gamma $Ga(g,h)$ distribution. Using Bayes Theorem, the posterior density is \begin{align*} \pi(\theta|\xvec)&\propto\pi(\theta)\,f(\xvec|\theta) \\ &\propto\frac{h^g\,\theta^{g-1}e^{-h\theta}}{\Gamma(g)}\times \theta^n e^{-n\overline{x^2}\theta},\quad\quad\theta>0, \end{align*} that is, $\theta|\xvec\sim Ga(G=g+n,H=h+n\overline{x^2})$. We saw in Example 1.8 that to make a gamma distribution vague, we take $g\to 0$ and $h\to 0$. Therefore using a vague prior distribution gives the posterior distribution as $\theta|\xvec\sim Ga(n,n\overline{x^2})$. \item[(ii)] The asymptotic posterior distribution (as $n\to\infty$) is \begin{equation*} \theta|\xvec\sim N\left(\hat\theta,J(\hat\theta)^{-1}\right), \end{equation*} where \begin{equation*} J(\theta)=-\frac{\partial^2}{\partial\theta^2} \log f(\xvec|\theta)=\frac{n}{\theta^2}. \end{equation*} Now \begin{eqnarray*} \frac{\partial}{\partial\theta} \log f(\xvec|\theta)=0 \quad\quad&\Longrightarrow&\quad\quad \frac{n}{\hat\theta}-n\overline{x^2}=0 \\ &\Longrightarrow&\quad\quad \hat\theta=1\left/\,\overline{x^2}\right. \\ &\Longrightarrow&\quad\quad J(\hat\theta)^{-1}=\frac{1}{n}\left(\overline{x^2}\right)^2. \end{eqnarray*} Therefore, for large $n$, the posterior distribution for $\theta$ is \begin{equation*} \theta|\xvec\sim N\left(1\left/\,\overline{x^2}\right.,\,\frac{1}{n}\left(\overline{x^2}\right)^2\right) \quad\quad\text{approximately}. \end{equation*} Note that the posterior distribution is centred on the likelihood mode and its variance tends to zero as $n\to\infty$. \end{itemize} \subsection*{Question 4} \begin{itemize} \item [(a)] We have: \begin{tabular}{|c|c|c|}\hline \centering $y$ &Cumulative probability &Decision? \\ \hline 5 &$0.45 $ &Keep going \\ 6 &$0.45+0.17=0.62$ &Keep going \\ 4 &$0.62+0.13=0.75$ &Keep going \\ 3 &$0.75+0.10=0.85$ &Keep going \\ 7 &$0.85+0.05=0.90$ &Stop \\ \hline \end{tabular} So the 90\% prediction set is $\{5, 6, 4, 3, 7\}$,in other words $\{3, 7\}$. \item [(b)] We have \[ E(y|x=4) = 1\times 0.01 + 2\times 0.04 + \cdots + 10\times 0.01 = 4.96. \] We also have \[ E(y^2|x=4) = 1^2 \times 0.01 + 2^2 \times 0.04 + \cdots + 10^2 \times 0.01 = 26.7. \] Therefore \[ Var(y|x=4) = 26.7 - 4.96^2 = 2.0984, \] giving $SD(y|x=4) = \sqrt{2.0984} = 1.4486$. \end{itemize} Therefore, \begin{eqnarray*} f(Y|x=4) &=&\binom{10}{y}\frac{B(y+9,33-y)}{B(9,23)} \\ & & \\ &=& \binom{10}{y} \frac{\frac{(y+9-1)!(33-y-1)!}{(y+9+33-y-1)!}}{\frac{(2-1)!(23-1)!}{(2+23-1)!}} \pause \\ & & \\ &=&\binom{5}{y} \times \frac{(y+1)!(27-y)!\times 24!}{22! \times 29!}. \end{eqnarray*} \pause So, for example, \[ P(Y=0|x=1) = \binom{5}{0} \times \frac{1! \times 27! \times 24!}{22! \times 29!} = 0.680 \quad \text{(See Table 4.4)} \] } } \section*{Distributions for data} \subsection*{Binomial distribution} If $X|\theta\sim Bin(k,\theta)$ then it has probability function \[ f(x|\theta)=\begin{pmatrix}k\\x\end{pmatrix}\theta^x(1-\theta)^{k-x} \,,\quad x=0,1,\ldots,k, \] where $k$ is a positive integer and $0<\theta<1$. Also, $E(X)=k\theta$ and $Var(X)=k\theta(1-\theta)$. \subsection*{Exponential distribution} If $X|\lambda\sim Exp(\lambda)$ then it has density \[ f(x|\lambda)=\lambda e^{-\lambda x}\,,\quad x>0 \] where $\lambda>0$. Also, $E(X)=1/\lambda$ and $Var(X)=1/\lambda^2$. \subsection*{Gamma distribution} If $X|\alpha,\lambda\sim Ga(\alpha,\lambda)$ then it has density \[ f(x|\alpha,\lambda)=\frac{\lambda^\alpha x^{\alpha-1}e^{-\lambda x}}{\Gamma(\alpha)}\,, \quad x>0, \] where $\alpha>0$ and $\lambda>0$. Also, $E(X)=\alpha/\lambda$ and $Var(X)=\alpha/\lambda^2$. \subsection*{Inverse Gaussian distribution} If $X|\mu,\lambda\sim IG(\mu,\lambda)$ then it has density \[ f(x|\mu,\lambda)=\frac{\sqrt\lambda}{\sqrt{2\pi x^3}} \exp\left\{-\frac{\lambda(x-\mu)^2}{2\mu^2 x}\right\}\,, \quad x>0, \] where $\mu>0$ and $\lambda>0$. Also, $E(X)=\mu$ and $Var(X)=\mu^3/\lambda$. \subsection*{Laplace distribution} If $X|\mu,\sigma\sim La(\mu,\sigma^2)$ then it has density \[ f(x|\mu,\sigma)=\frac{1}{\sqrt{2}\,\sigma} \exp\left\{-\frac{\sqrt{2}|x-\mu|}{\sigma}\right\}\,, \quad x>0, \] where $\mu\in\mathbb{R}$ and $\sigma>0$. Also, $E(X)=\mu$ and $Var(X)=\sigma^2$. \subsection*{Log-normal distribution} If $X|\mu,\sigma\sim LN(\mu,\sigma^2)$ then it has density \[ f(x|\mu,\sigma)= \frac{1}{x\sigma\sqrt{2\pi}} \exp\left\{-\frac{1}{2\sigma^2}(\log x-\mu)^2\right\},\quad x>0 \] where $\mu\in\mathbb{R}$ and $\sigma>0$. Also, $E(X)=e^{\mu+\sigma^2/2}$, $Var(X)=(e^{\sigma^2}-1)e^{2\mu+\sigma^2}$. Further $E(\log X)=\mu$ and $Var(\log X)=\sigma^2$. \subsection*{Negative binomial distribution} If $X|\theta\sim NegBin(k,\theta)$ then it has probability function \[ f(x|\theta)=\binom{x-1}{k-1}\theta^k(1-\theta)^{x-k},\quad x=k,k+1,\ldots~, \] where $k$ is a positive integer and $0<\theta<1$. Also $E(X)=k/\theta$ and $Var(X)=k(1-\theta)/\theta^2$. \subsection*{Normal distribution} If $X|\mu,\tau\sim N(\mu,1/\tau)$ then it has density \[ f(x|\mu,\tau)=\left(\frac{\tau}{2\pi}\right)^{1/2} \exp\left\{-\frac{\tau}{2}(x-\mu)^2\right\},\quad x\in\mathbb{R} \] where $\mu\in\mathbb{R}$ and $\tau>0$. Also, $E(X)=\mu$ and $Var(X)=1/\tau$. The distribution has the following quantiles \begin{center} \begin{tabular}{c|ccccc} $ x $ & 1.2816 & 1.6449 & 1.9600 & 2.3263 & 2.5758 \\ \hline $ \Pr(X\theta \] where $\alpha>0$ and $\theta>0$. Also $E(X)=\alpha\theta/(\alpha-1)$, $\alpha>1$ and $Var(X)=\alpha\theta^2/\{(\alpha-1)^2(\alpha-2)\}$, $\alpha>2$. \subsection*{Poisson distribution} If $X|\theta\sim Po(\theta)$ then it has probability function \[ f(x|\theta)=\displaystyle \frac{\theta^xe^{-\theta}}{x!}\,, \quad x=0,1,\ldots, \] where $\theta>0$. Also, $E(X)=\theta$ and $Var(X)=\theta$. \subsection*{Rayleigh distribution} If $X|\theta\sim R(\theta)$ then it has density \[ f(x|\theta)=2x\theta e^{-\theta x^2},\quad x>0, \] where $\theta>0$. Also $E(X)=\sqrt{\pi/(4\theta)}$ and $Var(X)=(4-\pi)/(4\theta)$. \subsection*{Uniform distribution} If $X|\phi,\theta\sim U(\phi,\theta)$ then it has density \[ f(x|\phi,\theta)=\frac{1}{\theta-\phi}\,,\quad \phi0$. Also $I_0(\cdot)$ is the modified Bessel function of order zero. \subsection*{Weibull distribution} If $X|\beta,\lambda\sim Wei(\beta,\lambda)$ then it has density \[ f(x|\beta,\lambda)=\beta\lambda^\beta x^{\beta-1}e^{-\lambda^\beta x^\beta}\,, \quad x>0, \] where $\beta>0$ and $\lambda>0$. Also, $E(X)=\Gamma(1+1/\beta)/\lambda$ and $Var(X)=\Gamma(1+2/\beta)/\lambda^2-E(X)^2$. \newpage \section*{Distributions for prior beliefs} \subsection*{Beta distribution} If $\theta\sim Beta(g,h)$ then it has density \[ \pi(\theta)=\frac{\theta^{g-1}(1-\theta)^{h-1}}{B(g,h)}\,, \quad 0<\theta<1, \] where $g>0$ and $h>0$. Also, $E(\theta)=g/(g+h)$ and $Var(\theta)=gh/\{(g+h)^2(g+h+1)\}$. \subsection*{Exponential distribution} If $\theta\sim Exp(h)$ then it has density \[ \pi(\theta)=h e^{-h\theta}\,,\quad \theta>0, \] where $h>0$. Also, $E(\theta)=1/h$ and $Var(\theta)=1/h^2$. \subsection*{Gamma distribution} If $\theta\sim Ga(g,h)$ then it has density \[ \pi(\theta)=\frac{h^g \theta^{g-1}e^{-h\theta}}{\Gamma(g)}\,, \quad \theta>0, \] where $g>0$ and $h>0$. Also, $E(\theta)=g/h$ and $Var(\theta)=g/h^2$. \subsection*{Generalised $\boldsymbol{t}$ distribution} \label{page:gt} If $\mu\sim t_a(b,c)$ then it has density \[ \pi(\mu)=\frac{\Gamma\left(\frac{a+1}{2}\right)} {\sqrt{ac\pi}\,\Gamma\left(\frac{a}{2}\right)} \left\{1+\frac{(\mu-b)^2}{ac}\right\}^{-\frac{a+1}{2}},\quad \mu\in\mathbb{R}, \] where $b\in\mathbb{R}$, $a>0$ and $c>0$. Also, $E(\mu)=b$ and $Var(\mu)=ac/(a-2)$ if $a\geq 2$. \subsection*{Inverse Chi distribution} \label{page:invchi} If $\sigma\sim\textit{Inv-Chi}(a,b)$ then it has density \[ \pi(\sigma|a,b)=\frac{2b^a\sigma^{-2a-1}e^{-b/\sigma^2}}{\Gamma(a)},\quad \sigma>0, \] where $a>0$, $b>0$ and $\Gamma(a)$ is the gamma function. Also $E(\sigma)=\sqrt{b}\,\Gamma(a-1/2)/\Gamma(a)$ and $Var(\sigma)=b/(a-1)-E(\sigma)^2$ if $a>1$. The name of the distribution comes from the fact that $1/\sigma^2\sim Ga(a,b)\equiv \chi^2_{2a}/(2b)$. \subsection*{Log-normal distribution} If $\theta\sim LN(b,c^2)$ then it has density \[ \pi(\theta)=\frac{1}{\sqrt{2\pi}\,c\,\theta} \exp\left\{-\frac{1}{2c^2}(\log\theta-b)^2\right\},\quad \theta>0 \] where $b\in\mathbb{R}$ and $c>0$. Also, $E(\theta)=e^{b+c^2/2}$, $Var(\theta)=(e^{c^2}-1)e^{2b+c^2}$. Further $\log\theta\sim N(b,c^2)$ and so $E(\log\theta)=b$ and $Var(\log\theta)=c^2$. \subsection*{Normal distribution} If $\mu\sim N(b,1/d)$ then it has density \[ \pi(\mu)=\left(\frac{d}{2\pi}\right)^{1/2} \exp\left\{-\frac{d}{2}(\mu-b)^2\right\},\quad \mu\in\mathbb{R}, \] where $b\in\mathbb{R}$ and $c>0$. Also, $E(\mu)=b$ and $Var(\mu)=1/d$. \subsection*{Normal-gamma distribution} If $\mutau\sim NGa(b,c,g,h)$ then it has density \[ \pi(\mu,\tau)\propto\tau^{g-\frac{1}{2}}\exp\left\{-\frac{\tau}{2} \left[c(\mu-b)^2+2h\right]\right\},\quad \mu\in\mathbb{R}, ~\tau>0 \] where $b\in\mathbb{R}$ and $c,g,h>0$. Also, $\mu|\tau\sim N\left(b,\dfrac{1}{c\tau}\right)$, $\tau\sim Ga(g,h)$ and has marginal distribution $\mu\sim t_{2g}\left(b,\dfrac{h}{gc}\right)$. \subsection*{Uniform distribution} If $\theta\sim U(a,b)$ then it has density \[ \pi(\theta)=\frac{1}{b-a}\,,\quad a<\theta