Past exam of the mathematics course of the University of Cambridge / 2013 / iii / Paper 34 (Applied Bayesian Statistics)

Each observation has uniform distribution density ${\theta }^{-1}{1}_{\{0<y_j<\theta \}}$. Independence makes the likelihood function their product. For positive observations the support reduces to $\theta >M$, givingFor data outside the positive orthant the likelihood is zero. Changing the convention at the endpoint does not change a continuous posterior distribution.

Multiplying the likelihood function by the Pareto distribution prior givesThe integral of this kernel is $m^{-(\alpha +n)}/(\alpha +n)$, so the normalized posterior distribution isThus it is $\mathrm{Pareto}(\alpha +n,m)$. This proves uniform-Pareto conjugacy: applying Bayes theorem preserves the family of prior distributions.

The joint prior predictive distribution is the Bayesian model evidence, obtained by integrating over the shared parameter:It is zero otherwise. This uniform-Pareto model evidence is a joint density, not the product of separately marginalized observation densities: mixing over the common parameter induces dependence.

Exchangeability means that the joint prior distribution is invariant under breed relabelling:for every permutation $\pi$. The labels carry no prior information about maximum size. This is reasonable for comparable breeds before observing their measurements when no covariates or biological knowledge distinguish them.

Exchangeability does not imply independence. In a hierarchical Bayesian model, the parameters can be independent conditional on a shared hyperparameter but dependent after it is integrated out. Known systematic biological differences would call for a model of those differences, with exchangeability only for the remaining unexplained variation.

Use the intended hierarchical Bayesian model: the breed parameters are conditionally independent given $\alpha ,\beta$, and observations are independent given their breed parameters. Put $m_i=\max (\beta ,M_i)$. The uniform-Pareto model evidence factorizes over breeds:Identical marginal prior distributions alone would not determine this product; conditional independence is the additional assumption.

The Type II maximum likelihood estimator maximizes the Bayesian model evidence after integrating out the breed parameters:For positive sample sizes define $S={\sum }_i\log (m_i/\beta )\ge 0$. Up to an additive constant, the log-likelihood isDifferentiation gives the interior score equationAlso ${\ell }^{\prime \prime }(\alpha )=-I/{\alpha }^2+{\sum }_i(\alpha +n_i{)}^{-2}<0$. If $S>0$, the derivative decreases from infinity to $-S$, so a unique finite maximum exists. When all sample sizes equal $n$, the equation becomes $In/[\alpha (\alpha +n)]=S$, givingFor $S=0$ no finite maximum exists. Plugging this hyperparameter estimate into the breed posterior distributions is an Empirical Bayes method.

Here $m_i=\beta$ for every breed, so $S=0$ andThe Bayesian model evidence increases towards its supremum as $\alpha \to \infty$. ThusFor a Pareto distribution, $\mathbb{P}(\theta >\beta (1+\epsilon ))=(1+\epsilon {)}^{-\alpha }\to 0$. The limiting empirical prior, and each corresponding posterior distribution, collapses onto $\beta$.

The boundary result is coherent within the assumed model: the observations favour the smallest allowed upper limits. Nevertheless, reporting exact concentration on a prespecified bound is overconfident for finite data and ignores hyperparameter uncertainty. A proper hyperprior on $\alpha$, sensitivity analysis for $\beta$, or a scientifically justified restriction on concentration avoids treating this limit as certain biological knowledge.

Fit two hierarchical Bayesian models to the same observations, one with the bounded uniform distribution sampling density and the other with an exponential distribution density. Specify rate versus mean in the exponential model and assign appropriate, scientifically comparable proper prior distributions; that parameter is no longer a literal maximum size. Obtain posterior distributions, for example by Markov chain Monte Carlo, and check convergence.

At the same observational level in both models use the Bayesian deviance $D(\varphi )=-2\log p(\mathbf{y}\mid \varphi )$, excluding prior densities. Retain the same likelihood constants and consistently either condition on breed effects or integrate them out. Estimate $\overline{D}=\mathbb{E}[D(\varphi )\mid \mathbf{y}]$, evaluate $D(\overline{\varphi })$ at the posterior mean, and compute the effective parameter count in DIC:Smaller DIC favours the fit–complexity tradeoff. Supplement it with posterior predictive checks; the uniform model's parameter-dependent support and the hierarchical structure make the criterion a diagnostic rather than an automatic definitive decision.

For a regular scalar statistical model, the one-observation Fisher information isUnder the usual differentiation and interchange-of-integral conditions it also equals $-{\mathbb{E}}_{\theta }[{\partial }_{\theta }^2\log p(Y\mid \theta )]$. The Jeffreys prior isThis prior measure is invariant under smooth one-to-one reparameterization, but may be improper. For independent identically distributed observations the information multiplier $n$ changes only the prior's proportionality constant.

For the scale family, substitute $u=y/\sigma$ in the expectation:This depends on $f,g$ but not on $\sigma$. It holds whenever the expectation exists, including an extended expectation for nonnegative $g$; an undefined difference of two infinite integrals is not an expectation.

Assume the scale family has the differentiability needed for Fisher information, and the following constant is finite and positive. With ${\ell }_f(u)=\log f(u)$, the score function isPart (b) therefore giveswhere $C$ is independent of $\sigma$. The Jeffreys prior for a scale parameter is consequentlyThe reciprocal-of-a-reciprocal in the TeX is a transcription error; the PDF has ${\sigma }^{-1}$. The expected-Hessian information identity should not be applied indiscriminately to nonregular parameter-dependent supports.

For $\tau =c\sigma$ with $c>0$, the change of variables givesEquivalently,This is the scale-invariant prior as a measure. Because its integral over $(0,\infty )$ diverges, it is an improper prior, not a normalized probability distribution on that range.

The normalized log-uniform distribution here has density $1/(x\log 10)$. A leading digit $i$ corresponds to $i\le X<i+1$, soThese are exactly the Benford law probabilities. Continuous densities make endpoint conventions immaterial.

The claim needs a complete number of logarithmic decades. For the normalized log-uniform distribution on the specified range, $U={\log }_{10}X$ is uniform on $(a,b)$. Digit $i$ corresponds to $\{U\}\in [{\log }_{10}i,{\log }_{10}(i+1))$. If $b-a=m$ is a positive integer, the integral of a period-one indicator over $(a,a+m)$ is $m$ times its integral over one period. Thus Benford law from logarithmic uniformity givesThis proves the intended case of integers $a<b$, and even permits noninteger $a$ when the span is an integer.

For arbitrary real $a<b$, the exact formula instead isOnly finitely many terms are nonzero. For a counterexample take $a=0$, $b={\log }_{10}2$: then $1<X<2$, so the leading digit is always one. That is not Benford law. The unrestricted range in the PDF needs this qualification.

A predictive discrepancy statistic is a specified measurable function $T$ of a possible dataset, chosen to detect a feature relevant to the null hypothesis. Under the fully specified null density $p_0$, compare the observed $T(\mathbf{x})$ with its reference distribution, derived analytically or from replicated datasets ${\mathbf{X}}^{\mathrm{rep}}\sim p_0$.

For discrepancies whose large values indicate disagreement, useA lower-tail or two-sided discrepancy requires the corresponding comparison. The checking function is a statistic, not a density; no unknown parameter is fitted when the null is fully specified.

Assuming independent digits under Benford law, put $N={\sum }_i{y}_i$ and $p_i={\log }_{10}(1+1/i)$. The count vector has a multinomial distribution, with expected counts $E_i=N{p}_i$. Suitable predictive discrepancy statistics include the Pearson chi-squared statistic and multinomial deviance:The zero-count terms of $T_G$ have limiting value zero. A Monte Carlo method gives a direct null comparison even when expected counts are small.

For example, an original R implementation is:Each column returned by rmultinom is a replicated count vector. R recycles the nine expected counts down each column. The add-one ratio is a Monte Carlo test estimate. In WinBUGS, alternatively generate a replicated dmulti vector using fixed $p_i,N$, compute its discrepancy and monitor exceedance of the observed value. The null does not estimate unknown digit probabilities. Dependence or selection in the accounts would require an appropriate simulation model.

Changing currency or monetary units multiplies amounts by a common constant without changing their generating mechanism. This motivates scale invariance of decimal significands of a generic significand distribution.

Multiplication by $c$ adds ${\log }_{10}c$ to the logarithm, rotating its fractional part. A probability distribution on the unit circle invariant under every rotation must be uniform: equal-length arcs have equal mass, and subdivision fixes that mass to arc length. Uniform fractional logarithms then give Benford law. Thus unit or currency invariance motivates uniform logarithmic mantissas. Restricted ranges, prescribed thresholds and rounding can prevent a particular collection from having this invariance.

Yes, the pattern merits checking. Digits two and three are visibly more frequent than the Benford law predictions, whereas several large digits are scarce and nine is absent. Rough agreement for one and four does not remove this pattern.

Use the predictive discrepancy statistic comparison in part (h), together with examination of the accounts' selection and constraints. A discrepancy is evidence against this digit model, not by itself evidence sufficient to establish fabrication. No numerical test calculation is needed for this qualitative assessment.

The multinomial distribution has probability mass functionIt is zero outside this count simplex. The multinomial coefficient counts the individual category sequences giving the same aggregate counts; the probability vector satisfies ${\sum }_k{p}_{ik}=1$.

Set ${\eta }_{ik}={\beta }_k^T{\mathbf{x}}_i$, so ${\eta }_{i1}=0$. The multinomial logistic regression implies $p_{ik}=e^{{\eta }_{ik}}{p}_{i1}$. Normalizing givesMultiply the multinomial likelihoods for conditionally independent groups. The multinomial coefficients are constant in $\beta$, leavingThe reference constraint ${\beta }_1=0$ identifies the coefficients: a common shift of every category coefficient otherwise leaves all probabilities unchanged. The PDF places the full sum inside the numerator exponent; the TeX breaks that expression.

Assume conditional independence of the Poisson distributions, and write $u_i={\mu }_{i1}$, $S_i={\sum }_k{e}^{{\eta }_{ik}}$. The correct Poisson mass function has $e^{-\mu }$; the positive sign in the PDF's reminder is erroneous. The group likelihood isMultiply by the shape–rate Gamma distribution prior $b^a{u}_i^{a-1}{e}^{-b{u}_i}/\Gamma (a)$ and integrate. The gamma integral givesIndependent baseline priors give the product over groups. For fixed $a,b$, the leading factors do not depend on $\beta$, proving the requested Gamma-integrated baseline Poisson likelihood.

The Poisson trick introduces ${\lambda }_i=\log {\mu }_{i1}$. With an exactly flat log-rate prior, the baseline prior measure is $d{u}_i/u_i$. For $n_i>0$,Thus flat log-rate marginalization gives the multinomial likelihood:Alternatively, independent Poisson counts conditional on their total are exactly multinomial distributions with probabilities $e^{{\eta }_{ik}}/S_i$. Locally uniform coefficient priors make the coefficient posterior proportional to this likelihood within their flat region.

The printed BUGS prior is a proper normal distribution on ${\lambda }_i$, with variance $s^2=10^5$, since BUGS uses precision as its second argument. It is broad but not exactly flat. Substitution $v=S_i{u}_i$ shows that, apart from a coefficient-independent constant, the integrated likelihood is the multinomial kernel times the Gaussian log-rate correction to the Poisson trickIt depends on the coefficients through $S_i$. Dominated convergence gives $R_i\to 1$ as $s^2\to \infty$, so the finite-variance code yields an approximation to the multinomial posterior, not exact equality. Large log rates can make prior sensitivity relevant.

The cell-wise factors form a standard Poisson regression, convenient for BUGS. They avoid an explicitly constrained count vector and permit scalar log-concave updates for intercepts and coefficients; an exactly flat-log-rate implementation also gives gamma conditional baseline rates. This can be computationally efficient, although the actual log-normal baseline prior is not gamma-conjugate and introduces $I$ nuisance intercepts. Efficiency depends on the update scheme.

The first coefficient is the posterior mean log odds of choosing invertebrates rather than fish for smaller alligators at Lake Hancock. It corresponds tofor those invertebrate-to-fish odds. This is not the absolute probability of choosing the second category, because the other categories also enter the normalization.

The size coefficient adds $-1.524$ to the invertebrate-to-fish log odds on changing to the larger class, holding lake fixed. Its common-across-lakes odds ratio isabout a 78 percent reduction in the relative odds. Exponentiating a posterior mean log odds is a geometric summary, not the arithmetic posterior mean of the odds.

There are $(K-1)p=4\cdot 5=20$ free category coefficients plus $I=8$ group intercepts, giving 28 free parameters in the fitted Poisson model. The effective parameter count in DIC, $26.8$, is close to 28 and slightly smaller, consistent with some regularization or incomplete information. Comparing it with 20 would omit the nuisance intercepts. A direct conditional multinomial fit has 20 coefficients but a different observational likelihood, since it conditions on the totals.

A scoring rule in the reward convention gives $S(P,x)$ for announcing the probability distribution $P$ and then observing $x$. Under genuine belief $Q$, the expected score isA proper scoring rule satisfies $s(Q,Q)\ge s(P,Q)$ for every admissible $P,Q$. A strictly proper scoring rule additionally requires equality only when the distributions coincide:The expected scores must be well defined. A loss convention reverses the inequalities, but the question uses rewards.

Let $q$ be the genuine rain probability and $p$ the reported probability. The linear probability score has expected rewardThis is linear in $p$, so its optimal report and expected reward areAt $q=1/2$ every report scores $1/2$ on average. Truthful reporting scores $q^2+(1-q{)}^2$, strictly below the optimum for interior $q\ne 1/2$. Thus the rule is not proper: it rewards maximal confidence in the more likely outcome. The expected total over days is ${\sum }_t\max (q_t,1-q_t)$.

Let $q$ be the true density and $p$ the announced density. The advantage of truthful reporting for the logarithmic scoring rule isThe Kullback-Leibler divergence is nonnegative, with equality exactly when the densities agree almost everywhere. For instance, Jensen inequality gives ${\mathbb{E}}_q[\log (p/q)]\le \log {\mathbb{E}}_q[p/q]\le 0$; equality requires a constant likelihood ratio on the true support and no remaining mass outside it. ConsequentlyThis proves strict propriety when the expected log scores are well defined, for example with a finite true expected log density. A forecast assigning zero density to a set of positive true probability has score $-\infty$ there and cannot outperform the truth. The hint's strict $>0$ must be corrected to $\ge 0$ to allow its equality case.

Model $m$ gives the sequential forecast $f_{mt}(y_t)=p_m(y_t\mid y_1,\dots ,y_{t-1})$. The probability chain rule yieldsApplication of Bayes theorem makes the factors successive posterior predictive distributions, and their product is the Bayesian model evidence. Under conditionally independent sampling it isTaking logarithms gives the prequential log score identityTherefore the Bayes factor isThis solves the second half of part (d), which is missing from the TeX. Proper priors and finite positive evidences are needed; unrelated improper-prior normalizing constants do not cancel.

Let $M\in \{1,2\}$ be one fixed model indicator with equal prior probabilities, and put $w_{m,t}=\mathbb{P}(M=m\mid y_1,\dots ,y_{t-1})$. Within each model the parameters are already updated using the same past data. The Bayesian model averaging forecast is $f_t(y)={\sum }_m{w}_{m,t}{f}_{mt}(y)$. Upon observing $y_t$, Bayes theorem givesTaking the ratio cancels the denominator and gives the prescribed update because $f_{mt}(y_t)=e^{L_{mt}}$. Iteration yieldsThe weights are posterior model probabilities and their mixture is the full Bayesian predictive density. The indicator is fixed across days, rather than choosing a fresh model independently each morning.

Integrate each supplied predictive density to get its cumulative distribution function $F_t$, then form the sequential probability integral transformFor a correct continuous one-step conditional predictive model, $\mathbb{P}(U_t\le u\mid y_1,\dots ,y_{t-1})=u$. Iterating this identity shows that the $U_t$ are independent uniform variables. A histogram or quantile plot checks uniformity; serial plots and autocorrelations check for temporal structure left unexplained by the forecasts.

Also inspect empirical coverage of central prediction intervals, tail exceedances and interval widths, assessing calibration together with sharpness. These checks need only the supplied forecasts and observations. Compare chosen predictive discrepancy statistics with simulated uniform reference sequences. A total log score alone is a relative reward and does not provide a universal absolute goodness-of-fit threshold.