Sequential Bayesian updating of mean and variance of normal distribution

Question

I am trying to write some code to learn the parameters of a normal distribution. I am new to this, and I have patched together the equations from various sources, which may be part of the problem. In particular, I seem to manage to update the prior with a batch of data, though I am having some difficulties interpreting the constant involved. However, I struggle to move from that to sequential updating, probably due to the same probelm with the constant.

Take a data vector $x = (x_i, ... , x_n)$, with $x_i \sim \mathcal{N}(\mu,\lambda^{-1})$, assuming exchangeable $x_i$, and where $\lambda \triangleq \frac{1}{\sigma^2}$ is the precision. I want to learn $p(\mu,\lambda | x) = p(x | \mu , \lambda)p(\mu,\lambda)$. Exploiting that $p(\mu,\lambda) = p(\mu | \lambda) p(\lambda)$, I can use the following normal-gamma prior: $$ p(\mu,\lambda) = \mathcal{N}(\mu|\mu_0,(n_0\lambda)^{-1}) Gam(\lambda|\alpha_0,\beta_0), $$ where $n_0$ should be some sort of scaling constant.

I should then be able to use the following equations to learn the model parameters (see here): $$ \alpha = \alpha_0 + \frac{n}{2}\\ \beta = \beta_0 + \frac{1}{2}\sum_{i=1}^{n}(x_i - \mu)^2 = \beta_0 + \frac{1}{2}\sum_{i=1}^{n}(x_i - \overline{x})^2 + \frac{1}{2}\frac{n n_0}{n + n_0}(\overline{x} - \mu_0)^2) \\ \lambda = \frac{\alpha}{\beta} \\ E[\mu|x,\lambda] =\frac{n\lambda}{n\lambda + n_0\lambda} \overline{x} + \frac{n_0\lambda}{n\lambda + n_0\lambda}\mu_0 = \mu_0 + \frac{n\lambda}{n\lambda + n_0\lambda} (\overline{x} - \mu_0), $$ where $\overline{x} = \frac{1}{n}\sum_i x_i$ is the sample mean. This seems to work on simulated data, but the results are very sensitive to my choice of $n_0$. The best results obtain for $n_0 = \frac{\alpha_0}{\beta_0} var(x)$. Hence my first question:

Q1: Is $n_0 = \frac{\alpha_0}{\beta_0} var(x)$ the correct way of fixing this `constant'? Can it be referred to as a constant at all, since it is clearly dependent on the data at hand?

Next, I would like to move on to sequential sampling, i.e. updating my prior with one observation at a time. The equation for $\alpha$ obtains trivially by setting $n=1$. For $\beta$, the middle term should drop out, since $x_i = \overline{x}$ for a single observation. My trouble is, once more, with $n_0$. For a single observation, the sample variance is not defined. This suggests that my solution above may have been wrong, or at least not generally applicable. I am also very confused of whether this `constant' should now be given one fixed value, or whether it should rather be adapted observation by observation.

Q2: how can I update the different parameters observation by observation? In particular, I am unclear what to do with the parameter $n_0$ if I set $n=1$.

Answer 1

The $n_0$ parameter is a parameter of a prior distribution. It is not estimated but assumed a priori. It serves as a "pseudocount" of the number of observations seen a priori and translates to "how strongly we believe" in the prior (Murphy, 2007). Notice this part of the equation

$$ \frac{n\lambda}{n\lambda + n_0\lambda} \overline{x} + \frac{n_0\lambda}{n\lambda + n_0\lambda}\mu_0 $$

where we take the weighted average between the observed $\overline{x}$ and the prior $\mu_0$. The weights depend on the number of samples you observed $n$ and $n_0$ "imaginary" samples from the prior. The more "imaginary" samples, the more weight goes to the prior.

If you want to update sequentially $n=1$ for each step. To choose the $n_0$ you need to consider how many samples would you consider to be enough to overwhelm the prior. Keep in mind that this is just intuition on how to think about it, because how strong the prior is depends on how much information it contains as compared to the data, so it is not only about the sample size.

Check "Conjugate Bayesian analysis of the Gaussian distribution" by Kevin P. Murphy for a more detailed description.

Answer 2

Thank you Tim for all the useful input and references!

Augmenting $n_0$, as I had suspected, is indeed what did the trick. With each observation, one additional "pseudocount" is thus added to the prior, making our belief in the prior increasingly strong. The notes of Kevin Murphy referenced by Tim above mention this case at the bottom of page 9. The equations thus become: $$ \alpha_i = \alpha_{i-1} + \frac{1}{2}\\ n_i = n_{i-1} + 1\\ \beta_i = \beta_{i-1} + \frac{1}{2}\frac{n_i}{n_i+1}(x_i - \mu_{i-1})^2 \\ \mu_i = \mu_{i-1} + \frac{1}{n_i + 1}(x_i - \mu_{i-1}) . $$