\documentclass[12pt,a4paper]{article} \usepackage{amsthm,amsmath,amsfonts} \usepackage{latexsym} \usepackage{fancyhdr} \usepackage{enumerate} %\usepackage{marginnote} \usepackage{graphicx} \usepackage{skull} \usepackage{hyperref} \usepackage{fancybox} %%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 1} \noindent Answers to the starred (\textbf{*}) questions should be submitted to Canvas by no later than \textcolor{red}{\bf 4pm, Friday 23rd February 2024}. Solutions to non-assessed questions will be made available in class and/or in Canvas. \subsection*{Question 1*} Nasa's new \textit{Earth Surface Mineral Dust Source Investigation} instrument investigates how airborne dust particles affect the climate. Since July 2020, the instrument has revealed that more than a billion metric tonnes -- the weight of 10,000 aircraft carriers -- of mineral dust are carried through the atmosphere from deserts and dry regions to other areas on earth. One area that has been particularly badly affected is southeastern Europe, in particular Greece, Cyprus and Turkey. In a joint project, data scientists at the Cyprus Institute and Newcastle University are interested in the probability $\theta$ that particulate matter in the air, during and immediatley following a dust storm, will result in an environmental health warning being issued for certain at-risk groups (i.e. people suffering from chronic lung conditions such as asthma). \begin{itemize} \item [(a)] After extensive discussions with a team of scientists -- including meteorologists and epidemiologists -- prior beliefs regarding $\theta$ were expressed according to the probability density function shown below in Figure 1. Look at it, and then answer the following questions. \begin{center} \begin{figure}[!h] \centering \includegraphics[width=9.5cm]{dust.pdf} \caption{Elicited prior density representing beliefs about $\theta$.} \end{figure} \end{center} \begin{itemize} \item [(i)] In plain English, briefly describe what Figure 1 tells us in terms of the project team's prior beliefs about $\theta$. \textcolor{gray}{\bf [2 marks]} \item [(ii)] The density shown in Figure 1 is that of a \textit{Beta distribution}. Choose the most appropriate prior from the list below, and clearly explain your choice: \textcolor{gray}{\bf [2 marks]} \[ \theta\sim Be(20, 77) \quad \text{\underline{or}} \quad \theta \sim Be(17, 5) \quad \text{\underline{or}} \quad \theta \sim Be(20, 20) \quad \text{\underline{or}} \quad \theta \sim Be(5, 17). \] \end{itemize} \item [(b)] Of the last ten dust storms to affect southeastern Europe, four resulted in environmental health warnings (in particular, the storm affecting Crete and Cyprus in summer 2021). \begin{itemize} \item [(i)] Assuming these storms occur independently, write down the likelihood function for $\theta$. \textcolor{gray}{\bf[2 marks]} \item [(ii)] Find the posterior distribution for $\theta$, and use plots/sketches and numerical summaries to help describe changes in beliefs about $\theta$ in light of this information. \textcolor{gray}{\bf[4 marks]} \end{itemize} \item [(c)]Find the posterior distribution for $\theta$ assuming an appropriate ignorance prior. How does this compare to the posterior obtained using the informative prior? \textcolor{gray}{\bf[3 marks]} \end{itemize} \subsection*{Question 2*} Suppose that a random sample $x_1,x_2,\ldots,x_n$ is obtained from a truncated unit exponential distribution, with density \[ f(x|\theta)=\begin{cases} e^{\theta-x},&x>\theta\\ 0,&\text{otherwise}, \end{cases} \] where $\theta\in\mathbb{R}$. Derive the posterior density for~$\theta$ assuming a normal $N(b,1/d)$ prior distribution. \textcolor{gray}{\bf [7 marks]} \\ Hint: \[ \theta^2-2b\theta - \frac{2n\theta}{d} \propto \left(\theta-\frac{db+n}{d}\right)^2. \] \subsection*{Question 3} Suppose you have a random sample $x_1, x_2, \ldots, x_n$. In the following models, derive the posterior density and name the posterior distribution (and its parameters). \begin{enumerate}[(i)] \item $f(x|\theta)=\theta^{x-1}(1-\theta )$, $x = 1,2,\ldots$ with a $Beta(3,2)$ prior distribution for~$\theta$. \item $f(x|\theta)=\displaystyle \frac{e^{-\theta}\theta^x}{x!}$, $x = 0,1,\ldots$ with a $Exp(2)$ prior distribution for $\theta$. \item $f(x|\theta)=\begin{pmatrix}4\\x\\ \end{pmatrix} \theta^x(1-\theta)^{4-x}$, $x=0,1,2,3,4$ with a $U(0,1)$ prior distribution for $\theta$. \end{enumerate} \subsection*{Question 4} Suppose that a random sample $x_1,x_2,\ldots,x_{10}$ is obtained from a uniform $U(0,\theta)$ distribution. Derive the posterior density for~$\theta$ assuming a gamma $Ga(20,1)$ prior distribution $\theta$. Hint: this distribution will depend on the maximum observed $x$-value $x_{\max}$. \subsection*{Question 5} The dimensions of a component from a long production run vary according to a $N(\mu ,1)$ distribution, and the mean dimension $\mu$ varies from production run to production run according to a $N(10,1/4)$ distribution. From one production run 12 components are drawn at random and their average dimension is found to be $10\frac{1}{3}$. On this information what is the probability that the mean component dimension is at least 10? \subsection*{Question 6} Suppose that the time in minutes required to serve a customer at a certain facility has an exponential distribution for which the parameter $\theta$ is unknown, and that the prior distribution of $\theta $ is a gamma distribution with mean 0.2 and standard deviation 1. If the average time required to serve a random sample of 20 customers is observed to be 3.8 minutes, determine the posterior distribution of $\theta$. \subsection*{Question 7} The data in Table 1 show the time intervals (in minutes) between eruptions of the ``Old Faithful'' geyser. Note that the average time between eruptions is 67.38 minutes, that is, 1.123 hours. If the eruptions occur randomly in time according to a Poisson process with a rate of $\theta$ eruptions per hour then the times between eruptions will form a random sample from an exponential distribution with rate~$\theta$. Making this assumption, establish a prior distribution for $\theta$ with $E(\theta) = 1$ and $Var(\theta) = 1$. Determine the posterior distribution for $\theta$. Comment on your analysis. \begin{table}[!h] \centering \begin{tabular}{|ccccccccccccccc|} \hline 70&64&72&76&80&48&88&53&71&56&69&72&76&54&76\\ 65&54&86&40&87&49&76&51&77&49&71&78&80&51&82\\ 49&80&43&83&49&75&47&78&71&69&63&64&82&68&71\\ 71&63&79&66&75&56&83&67&65&77&72&79&73&53&69\\ 53&78&55&67&68&73&53&70&69&66&79&48&90&49&78\\ 52&79&49&75&75&50&87&40&76&57&71&70&69&72&51\\ 84&43&73&73&70&84&71&79&58&73& & & & & \\ \hline \end{tabular} \caption{Time intervals between eruptions of the ``Old Faithful'' geyser.} \end{table} \newpage \section*{Solutions} \subsection*{Question 1} \begin{itemize} \item [(a)] \begin{itemize} \item [(i)] The prior mode is at $\theta = 0.2$, and most of the density is for $\theta<~0.6$, so it would appear that the epidemiologist thinks it unlikely that a dust storm will trigger an environmental health warning. \textcolor{red}{[2 marks for any similar sensible comment]} \item [(ii)] The correct prior is $\theta \sim Be(5, 17)$. Reasoning: \begin{itemize} \item The mode is clearly at 0.2. If $\theta \sim Be(a, b)$, then \[ Mode(\theta) = \frac{a-1}{a+b-2} = 0.2. \] For the four possibilities given, only $(a,b) = (20, 77)$ and $(a, b) = (5, 17)$ will satisfy this mode \item For these two possibilities, we have \[ \pi(\theta=0.2|a=20,b=77) = \frac{0.2^{19} \times 0.8^{76}}{B(20,77)} = \frac{0.2^{19} \times 0.8^{76}}{\frac{19!76!}{96!}} = 9.778, \] and \[ \pi(\theta=0.2|a=5,b=17) = \frac{0.2^4 \times 0.8^16}{B(5,17)} = \frac{0.2^4 \times 0.8^16}{\frac{4!16!}{21!}} = 4.582; \] comparing against the plot, we can see that $(a,b) = (5, 17)$ is clearly the correct choice. Of course, students could use R here to check: \begin{verbatim} > dbeta(0.2, 5, 17) [1] 4.582187 > dbeta(0.2, 20,77) [1] 9.778216 \end{verbatim} \end{itemize} \textcolor{red}{[2 marks for correct choice and reasoning; 0 marks otherwise -- even if they make the correct choice. Students could also simply plot all four options in R, and see which one matches that given in the question -- totally fine.]} \end{itemize} \item [(b)] \begin{itemize} \item [(i)] We assume $X|\theta \sim Bin(10, \theta)$, and so the likelihood term is \[ f(x=4|\theta) = \left(\begin{array}{c} 10 \\ 4 \end{array} \right) \theta^4 (1-\theta)^6, 0<\theta<1. \] \textcolor{red}{[2 marks; students must state the distribution they're assuming to get full marks]} \item [(ii)] We have \begin{eqnarray*} \pi(\theta|x=4) &\propto& \pi(\theta) \times f(x=4|\theta) \\ &\propto& \theta^{4}(1-\theta)^{16} \times \theta^4 (1-\theta)^6 \\ &\propto& \theta^{8}(1-\theta)^{22}, \quad 0<\theta<1, \end{eqnarray*} and so $\theta|x=4 \sim Be(9, 23)$. \textcolor{red}{[2 marks]} \newpage Students then need a couple of discussion points around the following plots/sketches and numerical summaries: \begin{center} \includegraphics[width=8cm]{ass1plot.pdf} \begin{tabular}{|c|c|c|c|} \hline $Be(a,b)$ &Prior, $(a=5, b=17)$ &Likelihood &Posterior, $(a=9, b=23)$ \\ \hline $Mode =\frac{a-1}{a+b-2}$ &0.2 &0.4 &0.267 \\ & & & \\ $SD =\sqrt{\frac{ab}{(a+b)^{2}(a+b+1)}}$ &0.087 &--- &0.078 \\ \hline \end{tabular} \end{center} \begin{itemize} \item In light of the data, we think it is more likely that a dust storm will trigger a health warning -- we can see this by an increase in the mode from prior to posterior \item Our uncertainty in $\theta$ has reduced slightly \end{itemize} \textcolor{red}{[2 marks -- one for each of the comments above. However, no marks should be awarded if the comments are not supported by plots/sketches and numerical summaries. No need to deduct marks if the table isn't exactly the same as the one above -- for example, students might not include the likelihood mode, which is fine.]} \end{itemize} \item [(c)] If we are prior ignorant, then we might just use $\theta \sim U(0,1)$, giving \begin{eqnarray*} \pi(\theta|x=4) &\propto & 1 \times \theta^4(1-\theta)^{6}, \end{eqnarray*} and so $\theta|x=4 \sim Be(5, 7)$. \textcolor{red}{[1 mark; fine if students use a Jeffreys prior instead]} This posterior is informed by the data only, and so (for example) \[ Mode(\theta|x=4) = \frac{4}{10} = 0.4 = \text{Likelihood mode}, \] compared to the informed posterior used earlier, whose mode (0.267) has been influenced heavily by the prior mode (0.2). \textcolor{red}{[2 marks for any sensible comment]} \end{itemize} \newpage \subsection*{Question 2} First, the likelihood: \begin{eqnarray*} f(\boldsymbol{x}|\theta) &=& \prod_{i=1}^{n} e^{\theta-x_{i}} \\ &=& e^{\theta+\theta \ldots +\theta} \times e^{-\sum_{i=1}^{n}x_{i}} \\ &=&e^{n\theta} \times e^{-n\bar{x}} \quad \text{\textcolor{red}{[2 marks for likelihood]}}\\ &\propto& e^{n\theta}, x_{\text{min}}>\theta. \end{eqnarray*} The prior is \[ \pi(\theta) = \left(\frac{d}{2\pi}\right)^{\frac{1}{2}} \exp\left\{-\frac{d}{2}\left(\theta-b\right)^{2}\right\}, \] giving \begin{eqnarray*} \pi(\theta|\boldsymbol{x}) &\propto& \exp\left\{n\theta -\frac{d}{2}(\theta-b)^2\right\} \\ &\propto& \exp\left\{n\theta -\frac{d}{2}(\theta^2 - 2b\theta +b^2)\right\} \\ &\propto& \exp\left\{n\theta-\frac{d}{2}\theta^2 +\frac{d}{2}2b\theta -\frac{d}{2}b^2\right\} \\ &\propto& \exp\left\{n\theta -\frac{d}{2}\theta^2+db\theta\right\}\times \exp\left\{-\frac{d}{2}b^2\right\} \\ &\propto& \exp\left\{n\theta -\frac{d}{2}\theta^2+db\theta\right\} \\ &\propto& \exp\left\{-\frac{d}{2}\left(\theta^2-\frac{2n\theta}{d}-2b\theta\right)\right\}. \end{eqnarray*} Using the hint, we see that \[ \pi(\theta|\boldsymbol{x}) \propto \exp\left\{-\frac{d}{2}\left(\theta-\frac{db+n}{d}\right)^{2}\right\}, \quad \theta0. \intertext{Therefore, using Bayes Theorem, the posterior density is} \pi(\theta|\xvec)&\propto\pi(\theta)\,f(\xvec|\theta) \\ &\propto e^{-2\theta}\times ke^{-n\theta}\theta^{n\bar x},\quad \theta>0 \\ &\propto\theta^{n\bar x}e^{-(n+2)\theta},\quad \theta>0, \end{align*} that is, $\theta|\xvec\sim Ga(n\bar x+1,n+2)$. \item [(iii)] The likelihood function is \[ f(\xvec|\theta)=\prod_{i=1}^n\begin{pmatrix} 4 \\ x_i \\\end{pmatrix} \theta^{x_i}(1-\theta)^{4-x_i} =k\theta^{n\bar x}(1-\theta)^{n(4-\bar x)} \] where $k$ is constant w.r.t. $\theta$. The prior distribution is $\theta\sim U(0,1)$, with density \[ \pi(\theta)=1,\quad 0<\theta<1. \] Therefore, using Bayes Theorem, the posterior density is \begin{align*} \pi(\theta|\xvec)&\propto\pi(\theta)\,f(\xvec|\theta) \\ &\propto 1\times k\theta^{n\bar x}(1-\theta)^{n(4-\bar x)}, \quad 0<\theta<1 \\ &\propto\theta^{n\bar x}(1-\theta)^{n(4-\bar x)},\quad 0<\theta<1, \end{align*} that is, $\theta|\xvec\sim Beta\{n\bar x+1,n(4-\bar x)+1\}$. \end{itemize} \subsection*{Question 4} The likelihood function is \begin{align*} f(\xvec|\theta)&=\prod_{i=1}^{10}\frac{1}{\theta},\quad x_1<\theta,x_2<\theta,\ldots,x_{10}<\theta\\ &=\theta^{-10},\quad x_{\max}<\theta. \end{align*} The prior distribution is $\theta\sim Ga(20,1)$, with density \begin{align*} \pi(\theta)&=\dfrac{1^{20}\theta^{20-1}e^{-\theta}}{\Gamma(20)},\quad \theta>0\\ &\propto\theta^{19}e^{-\theta},\quad \theta>0. \end{align*} Therefore, using Bayes Theorem, the posterior density is \begin{align*} \pi(\theta|\xvec)&\propto\pi(\theta)\,f(\xvec|\theta) \\ &\propto \theta^{19}e^{-\theta}\times\theta^{-10},\quad \theta>x_{\max}\\ &\propto \theta^9e^{-\theta},\quad \theta>x_{\max}. \end{align*} Note that this posterior distribution is a gamma $Ga(10,1)$ distribution, truncated below so that only values larger than $x_{\max}$ are possible. \subsection*{Question 5} From Example 1.3 (with $\tau=1$, $b=10$, $d=4$, $n=12$, $\bar x=31/3$), 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{40+124}{4+12}=10.25 \intertext{and} D&=d+n\tau=4+12=16=\frac{1}{0.25^2}, \end{align*} that is, $\mu|\xvec\sim N(10.25,0.25^2)$. Therefore $Pr(\mu>10|\xvec)=0.8413$ using\\ \texttt{1-pnorm(10,10.25,0.25)} \subsection*{Question 6} The gamma prior distribution for $\theta$ has mean 0.2 and standard deviation 1, whence $\theta\sim Ga(0.04,0.2)$. Also, we have $n=20$ observations with mean $\bar x=3.8$ from an $Exp(\theta)$ distribution. The likelihood function is \[ f(\xvec|\theta) =\prod_{i=1}^{20} \theta e^{-\theta x_i} \propto\theta^{20}e^{-20\bar x\theta} =\theta^{20}e^{-76\theta}. \] Using Bayes Theorem, the posterior density is \begin{align*} \pi(\theta|\xvec) &\propto\pi(\theta)\,f(\xvec|\theta) \\ &\propto\theta^{-0.96}e^{-0.2\theta}\times\theta^{20}e^{-76\theta}, \quad \theta>0 \\ &\propto\theta^{19.04}e^{-76.2\theta},\quad \theta>0, \end{align*} that is, $\theta|\xvec\sim Ga(G=20.04,H=76.2)$. This posterior distribution has mean 0.263 and standard deviation 0.059 (both to $3\,d.p.$). Clearly, the data has changed our (posterior) beliefs about the mean service rate~$\theta$: we now believe that $\theta$ takes a larger value (the mean has increased) and we are much more certain about likely values of $\theta$ (the variance is greatly reduced). \subsection*{Question 7} From the data we have $n=100$ and $\bar x=1.123$\,hours. If the data are a random sample from an exponential distribution then the mathematics used to derive the posterior distribution are more straightforward if we take a gamma prior distribution for $\theta$. (Note we could use any continuous distribution on the positive real line. Other possibilities include a log-normal distribution, a Weibull distribution and a Pareto distribution.) A gamma $Ga(g,h)$ prior distribution has mean $g/h$ and variance $g/h^2$ and so to make $E(\theta)=1$ and $Var(\theta)=1$ we need to take $g=h=1$. The likelihood function is \[ f(\xvec|\theta) =\prod_{i=1}^{100} \theta e^{-\theta x_i} \propto\theta^{100}e^{-100\bar x\theta} =\theta^{100}e^{-112.3\,\theta}. \] Using Bayes Theorem, the posterior density is \begin{align*} \pi(\theta|\xvec) &\propto\pi(\theta)\,f(\xvec|\theta) \\ &\propto\theta^0 e^{-\theta}\times\theta^{100}e^{-112.3\,\theta}, \quad \theta>0 \\ &\propto\theta^{100}e^{-113.3\,\theta},\quad \theta>0, \end{align*} that is, $\theta|\xvec\sim Ga(G=101,H=113.3)$. This distribution has mean 0.8914 and standard deviation 0.0887. The posterior mean is a little lower than the prior mean. Also the posterior standard deviation is considerably smaller than the prior standard deviation and so the data have been very informative about~$\theta$. Thus the posterior distribution will not be sensitive to minor changes to the prior mean and variance, that is, the posterior parameters $G$ and $H$ will be fairly similar for different (but similar) values of $g$ and~$h$. Also, given the amount of information in the data, it is likely that any other choice of prior distribution (with profile similar to the gamma distribution) will result in a very similar posterior distribution. Incidentally, by the Central Limit Theorem, this gamma posterior distribution is well approximated with a normal distribution with the same mean and variance. \end{document} \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