In this blog, I share my solutions on Mackay 2003 exercises chapter 27 page 342. Disclaimer: I don’t guarantee the validity of each answer. All of the answers are based on my own. However, I am also open for corrections. If you spot any flaws, feel free to contact me via email: arijonamarshal@gmail.com or marshal.arijona01@ui.ac.id. All of the notations follow the text book. For further details, please refer to https://www.inference.org.uk/itprnn/book.pdf.

Problem 1

A photon counter is pointed at a remote star for one minute, in order to infer the rate of photons arriving at the counter per minute, $λ$ . Assuming the number of photons collected r has a Poisson distribution with mean $λ$

$p (r | λ) = \exp (λ) \frac{λ^{r}}{r!}$

and assuming the improper prior P(λ) = 1/λ, make Laplace approximations to the posterior distribution

(a) over $λ$

(b) over $\log λ$ . [Note the improper prior transforms to $p (\log λ)$ is constant.

Problem 1a

First, let us compute the unnormalized posterior distribution $p^{*} (λ | r)$ and its log respectively:

\begin{array}{rcl} p^{*} (λ | r) & = & p (r | λ) p (λ) \\ = & \exp (- λ) \frac{λ^{r - 1}}{r!} \end{array}

$- \log p^{*} (λ | r) = λ - (r - 1) \log λ + \log r!$

Recall that we need the mode of distribution in order to perform Laplace approximation. We can compute the mode of $p^{*} (λ | r)$ via maximum a posteriori (MAP). MAP requires us to find $λ^{MAP}$ such that:

$λ^{MAP} = min_{λ} - \log p (λ | r)$

we can find MAP by deriving $- \log p^{*} (λ | r)$ , setting it to zero, and finally solving the equation for $λ$ :

\begin{array}{rcl} \frac{d - \log p^{*} (λ | r)}{d λ} & = & 1 - \frac{r - 1}{λ} = 0 \\ λ^{MAP} & = & r - 1 \end{array}

Now we can obtain the second order derivative $c$ on $λ = λ^{MAP}$ :

\begin{array}{rcl} c & = & {| - \frac{d^{2} \log p^{*} (λ | r)}{d λ^{2}} |}_{λ = λ^{MAP}} \\ = & - 1 (\frac{- (r - 1)}{λ^{2}}) \\ = & \frac{1}{r - 1} \end{array}

We are ready to construct our approximate distribution. First, let’s construct the unnnormalized approximation $q^{*} (λ)$ :

\begin{array}{rcl} q^{*} (λ) & \equiv & p^{*} (λ = λ^{MAP} | r) \exp (- \frac{c}{2} (λ - λ^{MAP})^{2}) \\ = & p^{*} (λ = λ^{MAP} | r) \exp (- \frac{1}{2 (r - 1)} (λ - (r - 1))^{2}) \end{array}

Our normalization factor $Z_{q}$ can be written as:

\begin{array}{rcl} Z_{q} & = & p^{*} (λ = λ^{MAP} | r) \sqrt{\frac{2 π}{c}} \\ = & p^{*} (λ = λ^{MAP} | r) \sqrt{2 π (r - 1)} \end{array}

Therefore, we obtain $q (λ) = \frac{1}{\sqrt{2 π (r - 1)}} \exp (- \frac{1}{2 (r - 1)} (λ - (r - 1))^{2}) = N (r - 1, r - 1)$ .

Problem 1b

From the description, $λ$ is transformed through the function $u (λ) = \log λ$ . Subsequently, the density of $λ$ is transformed to $p (u) = p (x) | \frac{d λ}{d u} | = p (x) λ$ . Using this rule, we obtain our unnormalized posterior $p^{*} (u (λ) | r)$ as follows:

\begin{array}{rcl} p^{*} (u (λ) | r) & = & p (u (λ)) p (r | u (λ)) \\ (1) & = & \exp (- λ) \frac{λ^{r}}{r!} λ = \exp (- λ) \frac{λ^{r + 1}}{r!} \end{array}

Observe that our unnormalized transformed posterior is just the transforming likelihood since our transformed prior is just a constant. Now, taking the log of $p^{*} (u (λ) | r)$ gives us:

$\begin{matrix} (2) & - \log p^{*} (u (λ) | r) = λ - (r + 1) \log λ + \log r! \end{matrix}$

Now, let’s derive $u (λ)^{MAP}$ :

\begin{array}{rcl} \frac{d - \log p^{*} (u (λ) | r)}{d u (λ)} & = & \exp (\log λ) - (r + 1) = 0 \\ u (λ)^{MAP} & = & \log (r + 1) \end{array}

We then derive our second derivative $c$ on $u (λ) = u (λ)^{MAP}$ :

\begin{array}{rcl} c & = & {| - \frac{d^{2} \log p^{*} (u (λ) | r)}{d u (λ)^{2}} |}_{u (λ) = u (λ)^{MAP}} \\ = & \exp (\log (r + 1)) \end{array}

We are ready to construct our approximate distribution. First, let’s construct the unnnormalized approximation $q^{*} (u (λ))$ :

\begin{array}{rcl} q^{*} (u (λ)) & \equiv & p^{*} (u (λ) = u (λ)^{MAP} | r) \exp (- \frac{c}{2} (u (λ) - u (λ)^{MAP})^{2}) \\ = & p^{*} (u (λ) = u (λ)^{MAP} | r) \exp (- \frac{(r + 1)}{2} (u (λ) - \log (r + 1))^{2}) \end{array}

Our normalization factor $Z_{q}$ can be written as:

\begin{array}{rcl} Z_{q} & = & p^{*} (u (λ) = u (λ)^{MAP} | r) \sqrt{\frac{2 π}{c}} \\ = & p^{*} (u (λ) = u (λ)^{MAP} | r) \sqrt{\frac{2 π}{r + 1}} \end{array}

Therefore, we obtain $q (u (λ)) = \frac{1}{\sqrt{\frac{2 π}{(r + 1)}}} \exp (- \frac{(r + 1)}{2} (u (λ) - \log (r + 1))^{2}) = N (\log (r + 1), \frac{1}{r + 1})$ .

Problem 2

Use Laplace’s method to approximate the integral

$Z (u_{1}, u_{2}) = \int_{- \infty}^{\infty} f (a)^{u_{1}} (1 - f (a))^{u_{2}} d a$

where $f (a) = 1 / (1 + e - a)$ and $u 1$ , $u 2$ are positive. Check the accuracy of the approximation against the exact answer (23.29, p.316) for $(u 1, u 2) = (1 / 2, 1 / 2)$ and $(u 1, u 2) = (1, 1)$ . Measure the error $(l o g Z_{p} - l o g Z_{q})$ in bits

Answer:

Let us define $p^{*} (a) = f (a)^{u_{1}} (1 - f (a))^{u_{2}}$ . Therefore, we have $\log p^{*} (a) = u_{1} \log f (a) + u_{2} \log (1 - f (a))$ . Next task is to determine the mode $a^{0}$ of $p^{*} (a)$ . Setting the derivative of $- \log p^{*} (a)$ to 0 and solving for $a$ , we obtain:

\begin{array}{rcl} \frac{d - \log p^{*} (a)}{d a} & = & - (\frac{u_{1}}{f (a)} f^{'} (a) - \frac{u_{2}}{1 - f (a)} f^{'} (a) = 0) \\ = & - (u_{1} (1 - f (a)) - u_{2} f (a)) \\ a^{0} & = & - \log \frac{u_{2}}{u_{1}} \end{array}

$f^{'} (a)$ refers to $\frac{d f (a)}{d a}$ . The second row comes from the fact that $f^{'} (a) = f (a) (1 - f (a))$ . Next, let us determine the second derivation.

\begin{array}{rcl} \frac{d^{2} - \log p^{*} (a)}{d a^{2}} & = & (u_{1} + u_{2}) f^{'} (a) \\ = & (u_{1} + u_{2}) f (a) (1 - f (a)) \end{array}

Now, we aim to evaluate $\frac{d^{2} - \log p^{*} (a)}{d a^{2}}$ at $a = a^{0}$ . Note that we have $f (a^{0}) = \frac{u_{1}}{u_{1} + u_{2}}$ and $1 - f (a) = \frac{u_{2}}{u_{1} + u_{2}}$ . Therefore

\begin{array}{r} c = {| \frac{d^{2} - \log p^{*} (a)}{d a^{2}} |}_{a = a^{0}} = (u_{1} + u_{2}) \frac{u_{1}}{u_{1} + u_{2}} \frac{u_{2}}{u_{1} + u_{2}} = \frac{u_{1} u_{2}}{u_{1} + u_{2}} \end{array}

The normalizing constant can be approximated by:

\begin{array}{rcl} Z_{p} ≃ Z_{q} & = & p^{*} (a^{0}) \sqrt{\frac{2 π}{c}} \\ = & {(\frac{u_{1}}{u_{1} + u_{2}})}^{u_{1}} {(\frac{u_{2}}{u_{1} + u_{2}})}^{u_{2}} \sqrt{\frac{2 π (u_{1} + u_{2})}{u_{1} u_{2}}} \end{array}

It’s time for evaluation !! for $(u_{1}, u_{2}) = (1, 1)$ , we have:

$Z_{q} (u_{1} = 1, u_{2} = 1) = \frac{1}{2} \times \frac{1}{2} \times 2 \sqrt{π} = \frac{\sqrt{π}}{2}$

$Z_{p} (u_{1} = 1, u_{2} = 1) = \frac{Γ (1) Γ (1)}{Γ (1 + 1)} = 1$

with the error:

$\log \frac{Z_{p}}{Z_{q}} = \log \frac{2}{\sqrt{π}} = 0.15 bit$

For the error measurement, we use base 2 logarithm. In other cases we use natural number.

while for $(u_{1}, u_{2}) = (\frac{1}{2}, \frac{1}{2})$ , we have:

$Z_{q} (u_{1} = 1 / 2, u_{2} = 1 / 2) = {\frac{1}{2}}^{1 / 2} \times {\frac{1}{2}}^{1 / 2} \times \sqrt{8 π} = \sqrt{2 π}$

$Z_{p} (u_{1} = 1 / 2, u_{2} = 1 / 2) = \frac{Γ (1 / 2) Γ (1 / 2)}{Γ (1)} = {1.77}^{2} = 3.313$

with the error: $\log \frac{Z_{p}}{Z_{q}} = \log \frac{3.313}{2.5} = 0.12 bit$

Problem 3

Linear regression. $N$ datapoints ${x^{n}, t^{n}}_{n = 1}^{N}$ are generated by the experimenter choosing each $x^{n}$ , then the world delivering a noisy version of the linear function

\begin{array}{rcl} y (x) & = & w 0 + w 1 x \\ t^{n} & \sim & N (y (x^{n}), σ_{ν}^{2}) \end{array}

Assuming Gaussian priors on $w_{0}$ and $w_{1}$ , make the Laplace approximation to the posterior distribution of $w_{0}$ and $w_{1}$ (which is exact, in fact) and obtain the predictive distribution for the next datapoint $t^{N + 1}$ , given $x^{N + 1}$

Answer :

Suppose that $w = (w_{0}, w_{1})$ , then our prior will be as follows:

\begin{array}{rcl} p (w) & = & N (w; 0, I) \\ \propto & \exp (- \frac{1}{2} w^{T} w) \\ \log p (w) & \propto & - \frac{1}{2} w^{T} w \end{array}

Suppose that we have $(x, t) = {x^{n}, t^{n}}_{n = 1}^{N}$ training-data. Since our likelihood is a Gaussian which depends on $x, w_{0},$ and $w_{1}$ , we obtain the likelihood as follows:

\begin{array}{rcl} p (t | x, w_{0}, w_{1}) & \propto & \prod_{n = 1}^{N} \exp (- \frac{1}{2 σ_{ν}^{2}} (t^{n} - y (x^{n}))^{2}) \\ = & \exp (- \frac{1}{2 σ_{ν}^{2}} \sum_{n = 1}^{N} (t^{n} - y (x^{n}))^{2}) \\ \log p (t | x, w_{0}, w_{1}) & \propto & - \frac{1}{2 σ_{ν}^{2}} \sum_{n = 1}^{N} (t^{n} - y (x^{n}))^{2} \end{array}

Having prior and likelihood. We are ready to compute the log of unnormalized posterior. We only need to apply the Bayes theorem to do so.

\begin{array}{rcl} \log p (w | x, t) & \propto & \log p (t | x, w) + \log p (w) \\ = & - \frac{1}{2 σ_{ν}^{2}} \sum_{n = 1}^{N} (t^{n} - y (x^{n}))^{2} - \frac{1}{2} w^{T} w \end{array}

Let’s define the unnormalized posterior as $p^{*} (w | x, t)$ . The rest of the solution is solving the Laplace approximation. The first step is to obtain $w^{MAP}$ via first derivation of $- \log p^{*} (w | x, t)$

\begin{array}{rcl} \frac{d - \log p^{*} (w | x, t)}{d w} & = & [\begin{array}{c} \frac{d - \log p^{*} (w | x, t)}{d w_{0}} \\ \frac{d - \log p^{*} (w | x, t)}{d w_{1}} \end{array}] \\ = & [\begin{array}{c} \frac{- 1}{σ_{ν}^{2}} \sum_{n = 1}^{N} [t^{n} - (w_{0} + w_{1} x^{n})] + w_{0} \\ \frac{- 1}{σ_{ν}^{2}} \sum_{n = 1}^{N} [t^{n} - (w_{0} + w_{1} x^{n})] x^{n} + w_{1} \end{array}] \end{array}

\begin{array}{rcl} w^{MAP} & = & [\begin{array}{c} w_{0}^{MAP} \\ w_{1}^{MAP} \end{array}] \\ = & \frac{1}{(n + σ_{ν}^{2}) (\sum_{n = 1}^{N} x^{n} + σ_{ν}^{2}) - (\sum_{n = 1}^{N} x^{n})^{2}} [\begin{array}{c} (\sum_{n = 1}^{N} (x^{n})^{2} + σ_{ν}^{2}) \sum_{n = 1}^{N} t^{n} - \sum_{n = 1}^{N} x^{n} \sum_{n = 1}^{N} x^{n} t^{n} \\ \sum_{n = 1}^{N} x^{n} \sum_{n = 1}^{N} t^{n} + (n + σ_{ν}^{2}) \sum_{n = 1}^{N} x^{n} t^{n} \end{array}] \end{array}

Next step is to obtain the Hessian matrix $H$ of $- \log p^{*} (w | x, t)$ :

\begin{array}{rcl} H & = & [\begin{array}{c} \frac{d^{2} - \log p^{*} (w | x, t)}{d w_{0}^{2}} & \frac{d^{2} - \log p^{*} (w | x, t)}{d w_{0} d w_{1}} \\ \frac{d^{2} - \log p^{*} (w | x, t)}{d w_{1} d w_{0}} & \frac{d^{2} - \log p^{*} (w | x, t)}{d w_{1}^{2}} \end{array}] \\ = & [\begin{array}{c} \frac{n}{σ_{ν}^{2}} + 1 & \frac{\sum_{n = 1}^{N} x^{n}}{σ_{ν}^{2}} \\ \frac{\sum_{n = 1}^{N} x^{n}}{σ_{ν}^{2}} & \frac{\sum_{n = 1}^{N} (x^{n})^{2}}{σ_{ν}^{2}} + 1 \end{array}] \end{array}

Now, we can approximate $p^{*} (w | x, t)$ with a distribution $q^{*} (w)$ : $\begin{array}{r} q^{*} (w) = p^{*} (w^{MAP}) \exp [\frac{- 1}{2} (w - w^{MAP})^{T} H (w - w^{MAP})] \end{array}$

Subsequently, we obtain the normalization factor $Z_{q}$

\begin{array}{r} Z_{q} = p^{*} (w^{MAP}) \sqrt{\frac{(2 π)^{K}}{det H}} \end{array}

with $K$ denotes the dimensionality of $w$ . Finally our approximate distribution follows a normal distribution

\begin{array}{rcl} q (w) & = & \frac{q^{*} (w)}{Z_{q}} \\ = & \frac{1}{\sqrt{\frac{(2 π)^{2}}{det H}}} \exp [\frac{- 1}{2} (w - w^{MAP})^{T} H (w - w^{MAP})] \\ = & N (w; w^{MAP}, H) \end{array}

Given a new data point $(x^{N + 1})$ , we compute the approximate predictive distribution:

\begin{array}{r} p (t^{N + 1} | x^{N + 1}) = \int p (t^{N + 1} | x^{N + 1}, w) q (w) d w \end{array}

However, compute this integral is often intractable. We can utilize Monte Carlo estimation to simplify the computation.

\begin{array}{r} p (t^{N + 1} | x^{N + 1}) \approx \sum_{i = 1}^{M} p (t^{N + 1} | x^{N + 1}, w^{i}), w^{i} \sim q (w) \end{array}

with M is the number of samples $w$ . The more samples we use the more accurate is the prediction.