# Normal distribution with unknown mean and variance

I'm reading about Bayesian Data Analysis by Gelman et al. and I came across with an example which starts like this:

Example. Normal distribution with unknown mean and variance.

We illustrate the approximate normal distribution with a simple theoretical example. Let $y_1, . . . , y_n$ be independent observations from a $N(\mu, \sigma^2)$ distribution, and, for simplicity, we assume a uniform prior density for $(\mu, \log \sigma)$. We set up a normal approximation to the posterior distribution of $(\mu, \log \sigma)$, which has the virtue of restricting $\sigma$ to positive values. To construct the approximation, we need the second derivatives of the log posterior density,

$$\log p(\mu, \log\sigma \,|\, y) = constant − n \log \sigma − \frac{1}{2\sigma^2} \left((n − 1)s^2 + n(\overline{y} − μ)^2\right),$$

where $\overline{y}=\frac1n\sum_{i=1}^n y_i$ and $s^2=\frac{1}{n-1}\sum_{i=1}^n(y_i-\overline{y})^2$.

This log of the posterior is escaping me for some reason...I'm not fully following how the author arrived in the expression above. I did try to calculate the posterior as the product of likelihood and prior but I didn't get the same result.

Is there a change of variables etc. happening under the hood in this example or what?

In summary my question is: How do we arrive to the log posterior expression? Steps?

For this you start by writing out your loglikelihood \begin{align} \log p(y|\mu,\sigma) &= \sum_i \log (y_i | \mu,\sigma) \\ &= \sum_i \left(-\frac{(y_i-\mu)^2}{2\sigma^2} - \log(\sqrt{2\pi\sigma^2}) \right)\\ &= -n\log \sigma - \frac{1}{2\sigma^2}\sum_i(y_i-\mu)^2 + \mbox{constant}. \end{align} Now what happens next is the following decomposition of the sum of squares term \begin{align} \sum_i (y_i - \mu)^2 &= \sum_i \left( y_i - \bar{y} + \bar{y} - \mu)\right)^2 \\ &= \sum_i \left((y_i-\bar{y})-(\mu-\bar{y})\right)^2, \end{align} you might like to go from here?

Ok if you are happy with that rearrangement of the likelihood then looking at the Bayesian aspects we are looking for the posterior $$p(\mu,\log \sigma | y) \propto p(y|\mu,\log \sigma) \pi(\mu,\log \sigma)$$ where $\pi(\cdot)$ is our prior, now this is assumed to be uniform which will also sometimes be given like $$\pi(\log \sigma) \propto 1,$$ now in specifying a uniform prior there is a little bit of handwaving going on - this is an improper prior if I want $\log \sigma$ to have an infinite support - I'm sure Gelman discusses this more in the book, a few releavent points are that improper priors can still lead to proper posteriors and we can compactify the infinite support if we want. Anyway assuming the $\mu$ also has an (independently of $\log\sigma$) uniform prior then taking the log of the posterior I have \begin{align} \log p(\mu,\log \sigma | y) &\propto \log p(y | \mu,\log \sigma ) + \log \pi(\log \sigma) + \log \pi(\mu) \\ &\propto \log p(y|\mu,\log \sigma) \end{align} Where to get the second line we used that both priors are proportional to a constant. Now as I mentioned as it stands this $\log \sigma$ on the righthand side is just a reparameterisation of the usual log-likelihood term which as far as writing this term out makes no difference, where it *will* make a difference is when we proceed to carry out differentiation to construct an approximation and so we will be differentiating with respect to $\log \sigma$ and not $\sigma$.

Finally why carry out this change of variables? Well I suspect what the author is going to do is construct a Laplacian approximation which amounts to a Gaussian approximation around the posterior mode, now if we were to construct this in the original $(\mu,\sigma)$ then the Guassian has support $\mathbb{R}^2$ and so will give positive probability to negative variances.

• Hi @Nadiels and thank you for your help! That already helped and I think I can finish yes after this. But until that, I have one more problem in your answer. In your first line you have $\log p(y|\mu,\sigma)$ but it shouldn't it be $\log p(y|\mu,\log \sigma)$? Because the log posterior is also expressed in terms of $\log \sigma$ and not just $\sigma$. This change of variables confuses me a bit x) – jjepsuomi Dec 7 '16 at 11:31
• Yeah I was just going to apply that change of variables at the end because this is just the most standard way of presenting the normal distribution, but $p(y|,mu,\sigma)$ is just as valid as $p(y|\mu,\log \sigma)$, where this change of variable would have an effect is in the prior $\pi(\sigma)$ compared to $\pi(\log \sigma)$, but in this case that is just left uniform – Nadiels Dec 7 '16 at 11:36
• Okay, that made it more clear. So you solve the loglikelihood and then you apply the change of variables which results in the form given by the author? If it's not a big effort for you, could you perhaps show just a little bit more how you would do this? This is a bit new subject for me and examples would be very nice to illustrate this. In any case, thank you very much for your valuable help! =) – jjepsuomi Dec 7 '16 at 11:47
• If I still haven't quite covered what you are interested in do say. When we reparameterise the log-likelihood term think of it at that stage as just another way of writing it, it does not affect the distribution of $y$, by which I mean $p(y|\sigma) = p(y|\log \sigma)$, but the importance of taking that transformation is it changes the support from $\sigma \in [0,\infty)$ to $\log \sigma \in \mathbb{R}$, which means that when we construct an approxmation $q(\log \sigma | y) \approx p(\log \sigma | y)$ and we say $q$ is Gaussian then the support of $q$ agrees with the support of $\log \sigma$ – Nadiels Dec 7 '16 at 12:40
• Hi @Nadiels, excellent answer, appreciate your effort! :) Sorry for taking time to answer, I had to be away from computer for a while. I think I got it now finally. In case I have more questions I will get back to the subject :) Thank you once again for your help! Superb! =) – jjepsuomi Dec 7 '16 at 14:59

Based on the answer by @Nadiels I'm trying to go through the math to arrive in the same log posterior as the author did in my question. I will continue from the answer given by @Nadiels:

So we have arrived at:

$$\log p(y|\mu,\sigma^2)= \text{constant}-n\log \sigma-\frac{1}{2\sigma^2}\sum_{i=1}^n ((y_i-\overline{y})-(\mu-\overline{y}))^2$$

$$=\text{constant}-n\log \sigma-\frac{1}{2\sigma^2}\sum_{i=1}^n \color{red}{ (y_i-\overline{y})^2}-\color{green}{2(y_i-\overline{y})(\mu-\overline{y})}+\color{blue}{(\mu-\overline{y})^2}$$

$$=\text{constant}-n\log \sigma-\frac{1}{2\sigma^2}\left(\color{red}{\sum_{i=1}^n (y_i-\overline{y})^2}-\color{green}{2(\mu-\overline{y})\sum_{i=1}^n(y_i-\overline{y})}+ \color{blue}{n(\overline{y}-\mu)^2} \right)$$

$$=\text{constant}-n\log \sigma-\frac{1}{2\sigma^2}\left(\color{red}{(n-1)s^2}-\color{green}{2(\mu-\overline{y})(n\overline{y}-n\overline{y})}+ \color{blue}{n(\overline{y}-\mu)^2} \right)$$

$$=\text{constant}-n\log \sigma-\frac{1}{2\sigma^2}\left((n-1)s^2+ n(\overline{y}-\mu)^2 \right),$$

so we got the same result as Gelman did in the book. Now this result indicates that we should have $p(y|\mu,\sigma^2)=p(y|\mu, \log \sigma)$. However, this is somewhat unclear to me since:

$$p(y|\mu,\sigma^2)=\prod_{i=1}^n\frac{1}{\sqrt{2\pi\sigma^2}}\exp\left(-\frac{1}{2\sigma^2}(y_i-\mu)^2\right)=\frac{1}{(2\pi\sigma^2)^{n/2}}\exp\left(-\frac{1}{2\sigma^2}\sum_{i=1}^n(y_i-\mu)^2\right),$$

and

$$p(y|\mu,\log\sigma)=\prod_{i=1}^n\frac{1}{\sqrt{2\pi\log\sigma}}\exp\left(-\frac{1}{2\log\sigma}(y_i-\mu)^2\right)$$

$$=\frac{1}{(2\pi\log\sigma)^{n/2}}\exp\left(-\frac{1}{2\log\sigma}\sum_{i=1}^n(y_i-\mu)^2\right),$$

So now to me it seems that $p(y|\mu,\sigma^2)\neq p(y|\mu,\log\sigma)$. Take for example $\sigma=e,$

$$p(y|\mu,e^2)=\frac{1}{(2\pi e^2)^{n/2}}\exp\left(-\frac{1}{2e^2}\sum_{i=1}^n(y_i-\mu)^2\right)$$

$$p(y|\mu,\log e)= p(y|\mu,1)=\frac{1}{(2\pi)^{n/2}}\exp\left(-\frac{1}{2}\sum_{i=1}^n(y_i-\mu)^2\right)$$

• Your confusion in the latter part is that when we consider some invertible parameter transform of a parameter $\theta$, say $\psi = g(\theta)$, then if we have a density $p(y|\theta)$ we can rewrite it as $p(y|\psi)$ by which we really mean $p(y|\theta = g^{-1}(\psi))$ - look at your density with $\sqrt{\log \sigma}$, $\log \sigma$ will be negative for $\sigma < 1$! – Nadiels Dec 8 '16 at 11:40
• So that taking $p(y| \mu, \log \sigma = 1)$ what we mean is that this is the density of Gaussian random variable $Y$ with mean $\mu$ and standard deviation $\sigma = \exp(\log \sigma)) = e$ – Nadiels Dec 8 '16 at 11:43
• I think I got it now. This notation confused me I guess, I'm used to the notation $p(y|\mu,x)$ where whatever $x$ is, it IS the scale parameter of the distribution of random variable $Y$. So $p(y|\mu,\log\sigma)$ was kinda like an indirect way to limit the values of $\sigma$. By the way, does the change of variables play a role in this, i.e. if $\phi=h(\theta)$ then $p(\phi)=p(\theta)\left|\frac{d\theta}{d\phi}\right|$? Is this connected to the subject we've been discussing? – jjepsuomi Dec 8 '16 at 11:54
• The change of variables does not play a role in this when it comes to the point wise evaluation of the conditional density of $y$, there are two important situation when it does make a difference however; 1) If we had specified our prior $\pi(\sigma)$ first and then decided to carry out $\log \sigma$ we would have to carry out the change of variables for the prior - that is avoided here by specifying the prior on $\log \sigma$. 2) If we wish to differentiate likelihood then we must differentiation with respect to our parameter $\log \sigma$, either directly or with the chain rule. – Nadiels Dec 8 '16 at 12:02
• You are very welcome – Nadiels Dec 8 '16 at 12:23