The Gaussian distribution

1.1 Moments form

The standard (univariate) Gaussian pdf is $$\psi :x\mapsto \frac{1}{\sqrt{2\pi }}\text{exp}(-\frac{x^2}{2}).$$ Typically we allow a scale-location parameterised version $$\varphi (x;\mu ,{\sigma }^2)=\frac{1}{\sqrt{2\pi {\sigma }^2}}{e}^{-\frac{(x-\mu {)}^2}{2{\sigma }^2}}$$ We call the CDF $$\mathrm{\Psi }:x\mapsto {\int }_{-\mathrm{\infty }}^x\psi (t)dt.$$ In the multivariate case, where the covariance $\mathrm{\Sigma }$ is strictly positive definite, we can write a density of the general normal distribution over ${\mathbb{R}}^k$ as $$\psi (x;\mu ,\mathrm{\Sigma })=(2\pi {)}^{-\frac{k}{2}}det(\mathrm{\Sigma }{)}^{-\frac{1}{2}}\exp (-\frac{1}{2}(x-\mu {)}^{\mathrm{\top }}{\mathrm{\Sigma }}^{-1}(x-\mu ))$$ If a random variable $Y$ has a Gaussian distribution with parameters $\mu ,\mathrm{\Sigma }$, we write $$Y\sim \mathcal{N}(\mu ,\mathrm{\Sigma })$$

1.2 Canonical form

Hereafter we assume that we are talking about multivariate Gaussians, or this is not very interesting. Davison and Ortiz (2019) recommends the “information parameterization” of Eustice, Singh, and Leonard (2006), which looks like the canonical form to me. $$p\left(\mathbf{x}\right)\propto e^{[-\frac{1}{2}{\mathbf{x}}^{\mathrm{\top }}\mathrm{\Lambda }\mathbf{x}+{\mathit{\eta }}^{\mathrm{\top }}\mathbf{x}]}$$ The parameters of the canonical form $\mathit{\eta },\mathrm{\Lambda }$ relate to the moments via $$\begin{array}{rl}\mathit{\eta }& =\mathrm{\Lambda }\mathit{\mu }\\ \mathrm{\Lambda }& ={\mathrm{\Sigma }}^{-1}\end{array}$$ The information form is convenient as multiplication of distributions (when we update) is handled simply by adding the information vectors and precision matrices — see Marginals and conditionals.

1.3 Tempered

Call the $\tau$-tempering of a density $p$ the density $p^{\tau }$ for $\tau \in \mathbb{R}$. As with all exponential families, it comes out simple if we use canonical form. For a multivariate Gaussian distribution in canonical (information) form, the density is expressed as

$$p(x)\propto \exp (-\frac{1}{2}{x}^{\mathrm{\top }}\mathrm{\Lambda }x+x^{\mathrm{\top }}\eta ),$$

where

• $\mathrm{\Lambda }$ is the precision matrix (the inverse of the covariance matrix $\mathrm{\Sigma }$, i.e., $\mathrm{\Lambda }={\mathrm{\Sigma }}^{-1}$),
• $\eta =\mathrm{\Lambda }\mu$ is the information vector,
• $\mu$ is the mean vector.

When we temper this density by $\tau$, we get

$$p_{\tau }(x)\propto {[\exp (-\frac{1}{2}{x}^{\mathrm{\top }}\mathrm{\Lambda }x+x^{\mathrm{\top }}\eta )]}^{\tau }=\exp (-\frac{1}{2}\tau x^{\mathrm{\top }}\mathrm{\Lambda }x+\tau x^{\mathrm{\top }}\eta ).$$

This expression is still in the form of a Gaussian distribution but with scaled parameters. Specifically, both the precision matrix and the information vector are scaled by $\tau$:

• New Precision Matrix: ${\mathrm{\Lambda }}^{\prime }=\tau \mathrm{\Lambda }$,
• New Information Vector: ${\eta }^{\prime }=\tau \eta$.

Thus, the tempered Gaussian retains the same structure but with parameters adjusted proportionally to $\tau$.

In the moments form, tempering a multivariate Gaussian distribution by a scalar $\tau \in (0,1]$ results in:

• Unchanged Mean: ${\mu }^{\prime }=\mu$
• Scaled Covariance: ${\mathrm{\Sigma }}^{\prime }={\displaystyle \frac{1}{\tau }}\mathrm{\Sigma }$

This means the tempered Gaussian has the same center but is more spread out, reflecting the reduction in “confidence” due to tempering.

2.1 What is Erf again?

This erf, or error function, is a rebranding and reparameterisation of the standard univariate normal cdf. It is popular in computer science, AFAICT because calling something an “error function” provides an amusing ambiguity to pair with the one we are used to with the “normal” distribution. There are scaling factors tacked on.

$$\mathrm{erf}(x)=\frac{1}{\sqrt{\pi }}{\int }_{-x}^x{e}^{-t^2}dt$$ which is to say $$\begin{array}{rl}\mathrm{\Phi }(x)& =\frac{1}{2}[1+\mathrm{erf}\left(\frac{x}{\sqrt{2}}\right)]\\ \mathrm{erf}(x)& =2\mathrm{\Phi }(\sqrt{2}x)-1\end{array}$$

2.2 Taylor expansion of $e^{-x^2/2}$

$$e^{-x^2/2}=\sum _{k=0}^{\mathrm{\infty }}(2^{(}-k)(-x^2{)}^k)/(k!).$$

2.3 Score

$$\begin{array}{rl}{\mathrm{\nabla }}_x\log \psi (x;\mu ,\mathrm{\Sigma })& ={\mathrm{\nabla }}_x(-\frac{1}{2}(x-\mu {)}^{\mathrm{\top }}{\mathrm{\Sigma }}^{-1}(x-\mu ))\\ & =-(x-\mu {)}^{\mathrm{\top }}{\mathrm{\Sigma }}^{-1}\end{array}$$

3.1 Stein’s lemma

Meckes (2009) explains Stein (1972)’s characterisation:

The normal distribution is the unique probability measure $\mu$ for which $$\int [f^{\prime }(x)-xf(x)]\mu (dx)=0$$ for all $f$ for which the left-hand side exists and is finite.

This is incredibly useful in probability approximation by Gaussians where it justifies Stein’s method.

3.2 ODE representation for the univariate density

$$\begin{array}{rl}{\sigma }^2{\varphi }^{\prime }(x)+\varphi (x)(x-\mu )& =0,\text{ i.e.}\\ L(x)& =({\sigma }^{2}D+x-\mu )\end{array}$$

With initial conditions

$$\begin{array}{rl}\varphi (0)& =\frac{e^{-{\mu }^2/(2{\sigma }^2)}}{\sqrt{2{\sigma }^2\pi }}\\ {\varphi }^{\prime }(0)& =0\end{array}$$

🏗 note where I learned this.

3.3 ODE representation for the univariate icdf

From (Steinbrecher and Shaw 2008) via Wikipedia.

Write $w:={\mathrm{\Psi }}^{-1}$ to keep notation clear. It holds that $$\frac{d^{2}w}{d{p}^2}=w{\left(\frac{dw}{dp}\right)}^2.$$ Using initial conditions $$\begin{array}{rl}w\left(1/2\right)& =0,\\ w^{\prime }\left(1/2\right)& =\sqrt{2\pi }.\end{array}$$ we generate the icdf.

3.4 Density PDE representation as a diffusion equation

Z. I. Botev, Grotowski, and Kroese (2010) notes

$$\begin{array}{rl}\frac{\partial }{\partial t}\varphi (x;t)& =\frac{1}{2}\frac{{\partial }^2}{\partial x^2}\varphi (x;t)\\ \varphi (x;0)& =\delta (x-\mu )\end{array}$$

Look, it’s the diffusion equation of Wiener process. Surprise! If you think about this for a while you end up discovering Feynman-Kac formulae.

3.5 Roughness of the density

Univariate —

$$\begin{array}{rl}{\Vert \frac{d}{dx}{\varphi }_{\sigma }\Vert }_2& =\frac{1}{4\sqrt{\pi }{\sigma }^3}\\ {\Vert {\left(\frac{d}{dx}\right)}^n{\varphi }_{\sigma }\Vert }_2& =\frac{\prod _{i<n}2n-1}{2^{n+1}\sqrt{\pi }{\sigma }^{2n+1}}\end{array}$$

4.1 Mills ratio

Mills’ ratio is $(1-\mathrm{\Phi }(x))/\varphi (x)$ and is a workhorse for tail inequalities for Gaussians. See the review and extensions of classic results in Dümbgen (2010), found via Mike Spivey. Check out his extended justification for the classic identity

$${\int }_x^{\mathrm{\infty }}\frac{1}{\sqrt{2\pi }}{e}^{-t^2/2}dt\le {\int }_x^{\mathrm{\infty }}\frac{t}{x}\frac{1}{\sqrt{2\pi }}{e}^{-t^2/2}dt=\frac{e^{-x^2/2}}{x\sqrt{2\pi }}.$$

7.1 Product of densities

A workhorse of Bayesian statistics is the product of densities, and it comes out in an occasionally useful form for Gaussians.

Let ${\mathcal{N}}_{\mathbf{x}}(\mathit{\mu },\mathbf{\Sigma })$ denote a density of $\mathbf{x}$, then $${\mathcal{N}}_{\mathbf{x}}({\mathit{\mu }}_1,{\mathbf{\Sigma }}_1)\cdot {\mathcal{N}}_{\mathbf{x}}({\mathit{\mu }}_2,{\mathbf{\Sigma }}_2)\propto {\mathcal{N}}_{\mathbf{x}}({\mathit{\mu }}_c,{\mathbf{\Sigma }}_c)$$ where $$\begin{array}{rl}{\mathbf{\Sigma }}_c& ={({\mathbf{\Sigma }}_1^{-1}+{\mathbf{\Sigma }}_2^{-1})}^{-1}\\ {\mathit{\mu }}_c& ={\mathbf{\Sigma }}_c({\mathbf{\Sigma }}_1^{-1}{\mathit{\mu }}_1+{\mathbf{\Sigma }}_2^{-1}{\mathit{\mu }}_2).\end{array}$$

Translating that to canonical form, with ${\mathcal{N}}_{\mathbf{x}}^{-1}(\mathit{\eta },\mathbf{\Sigma })$ denoting a canonical-form density of $\mathbf{x}$, then $${\mathcal{N}}_{\mathbf{x}}^{-1}({\mathit{\eta }}_1,{\mathbf{\Lambda }}_1)\cdot {\mathcal{N}}_{\mathbf{x}}^{-1}({\mathit{\eta }}_2,{\mathbf{\Lambda }}_2)\propto {\mathcal{N}}_{\mathbf{x}}^{-1}({\mathit{\eta }}_c,{\mathbf{\Lambda }}_c)$$ where $$\begin{array}{rl}{\mathbf{\Lambda }}_c& ={\mathbf{\Lambda }}_1+{\mathbf{\Lambda }}_2\\ {\mathit{\eta }}_c& ={\mathit{\eta }}_1+{\mathit{\eta }}_2.\end{array}$$

7.2 Orthogonality of conditionals

See pathwise GP.

13.1 In canonical form

Weirdly no-one seems to write out the canonical form of the conjugate prior for the Gaussian. For reasons of, presumably, tradition, weird ad hoc parameterisations such as mean-precision or mean variance are used. These are not convenient.

We start from a Gaussian likelihood with in canonical form. To put the Gaussian likelihood into exponential family form, let’s first consider the precision $ = 1/^2.$ The Gaussian distribution can be parameterised by the natural parameters $\eta =\left[\begin{array}{c}\rho \tau \end{array}\right]$, where $\rho =\mu \tau ,$ is the mean scaled by the precision, and $\tau$ is the precision.

Sanity check: Write out moment-form Gaussian and canonicalize. $$\begin{array}{rl}p(x\mid \mu ,{\sigma }^2)& =\frac{1}{\sqrt{2\pi {\sigma }^2}}{e}^{-\frac{(x-\mu {)}^2}{2{\sigma }^2}}\\ & =\sqrt{\frac{\tau }{2\pi }}{e}^{-\frac{1}{2}\tau (x-\rho /\tau {)}^2}\\ & =\sqrt{\frac{\tau }{2\pi }}{e}^{-\frac{1}{2}\tau (x^2-2\rho x/\tau +{\rho }^2/{\tau }^2)}\\ & =\frac{\sqrt{\tau }}{\sqrt{2\pi }}{e}^{-\frac{1}{2}\tau x^2+\rho x-\frac{1}{2}{\rho }^2/\tau }\\ & =\frac{1}{\sqrt{2\pi }}{e}^{-\frac{1}{2}\tau x^2+\rho x+\frac{1}{2}\log \tau -\frac{1}{2}{\rho }^2/\tau }\\ & =\frac{1}{\sqrt{2\pi }}\exp (\left[\begin{array}{c}x\\ -\frac{1}{2}{x}^2\end{array}\right]\cdot \left[\begin{array}{c}\rho \\ \tau \end{array}\right]+\frac{1}{2}\log \tau -\frac{1}{2}\frac{{\rho }^2}{\tau })\end{array}$$ The sufficient statistics $T(x)$ that correspond to these natural parameters are $T(x)=\left[\begin{array}{c}x-\frac{x^2}{2}\end{array}\right].$

The base measure $h(x)$ is $$h(x)=\frac{1}{\sqrt{2\pi }}.$$

The log-partition function $A(\eta )$ is $$A(\eta )=\frac{1}{2}\log \tau -\frac{1}{2}\frac{{\rho }^2}{\tau }.$$ Putting all these together, the canonical form of the Gaussian likelihood with respect to the natural parameters and sufficient statistics becomes: $$p(x\mid \rho ,\tau )=\frac{1}{\sqrt{2\pi }}\exp (\rho x-\frac{\tau x^2}{2}-(\frac{{\rho }^2}{2\tau }-\frac{1}{2}\log \tau ))$$

What is the conjugate prior for this? I am going to punt that making it Normal-Gamma will come out ok, why not? But I suspect we want the Normal term in the Normal-Gamma to be canonical form (?).

$$\begin{array}{rl}p\left(\left[\begin{array}{c}\rho \\ \tau \end{array}\right]\right|r,t,\kappa ,\gamma )& =\mathrm{NormalGamma}(\rho ,\tau |r,t,\kappa ,\gamma )\\ & =\mathcal{N}(\rho |r,t)\mathrm{Gamma}(\tau |\kappa ,\gamma )\\ & =\sqrt{\frac{t}{2\pi }}{e}^{-\frac{t}{2}({\rho }^2-2\rho r/t+r^2/t^2)}\frac{(-\gamma {)}^{\kappa +1}}{\mathrm{\Gamma }(\kappa +1)}{\tau }^{\kappa }{e}^{\tau \gamma }\\ & =\sqrt{\frac{t}{2\pi }}\frac{(-\gamma {)}^{\kappa +1}}{\mathrm{\Gamma }(\kappa +1)}{e}^{-\frac{1}{2}{r}^2/t}{\tau }^{\kappa }{e}^{-\frac{t}{2}{\rho }^2+\rho r+\tau \gamma }\\ & =\sqrt{\frac{t}{2\pi }}\frac{(-\gamma {)}^{\kappa +1}}{\mathrm{\Gamma }(\kappa +1)}{e}^{-\frac{1}{2}{r}^2/t}{e}^{\kappa \log \tau -\frac{t}{2}{\rho }^2+\rho r+\tau \gamma }\\ & =\left(\sqrt{\frac{t}{2\pi }}{e}^{-\frac{r^2}{2t}}\right)\left(\frac{(-\gamma {)}^{\kappa +1}}{\mathrm{\Gamma }(\kappa +1)}\right)\exp (\left[\begin{array}{c}\gamma \\ \kappa \\ r\\ t\end{array}\right]\cdot \left[\begin{array}{c}\tau \\ \log \tau \\ \rho \\ -\frac{{\rho }^2}{2}\end{array}\right])\end{array}$$ NB, $\gamma =-\beta$ in the classic Gamma rate param and $\kappa =k-1$. $\gamma$ is always negative.

Let us directly apply Bayes’ theorem. $$\begin{array}{rl}p(\rho ,\tau |x)& \propto p(x|\rho ,\tau )p(\rho ,\tau )\\ & =e^{\rho x-\frac{\tau x^2}{2}-(\frac{{\rho }^2}{2\tau }-\frac{1}{2}\log \tau )}\sqrt{\frac{t}{2\pi }}\frac{(-\gamma {)}^{\kappa +1}}{\mathrm{\Gamma }(\kappa +1)}{e}^{-\frac{1}{2}{r}^2/t}{e}^{\kappa \log \tau -\frac{t}{2}{\rho }^2+\rho r+\tau \gamma }\\ & \propto \sqrt{\frac{t}{2\pi }}\frac{(-\gamma {)}^{\kappa +1}}{\mathrm{\Gamma }(\kappa +1)}{e}^{\tau \gamma -\frac{\tau x^2}{2}+\frac{1}{2}\log \tau +\kappa \log \tau +\rho x+\rho r-\frac{t}{2}{\rho }^2-\frac{{\rho }^2}{2\tau }-\frac{1}{2}{r}^2/t}\\ & \propto e^{((\gamma -\frac{x^2}{2})\tau +(\frac{1}{2}+\kappa )\log \tau +(x+r)\rho x+(t+\frac{1}{\tau })(-\frac{{\rho }^2}{2})}\\ & =\mathrm{NormalGamma}(\rho ,\tau |r+x,t+\frac{1}{\tau },\kappa +\frac{1}{2},\gamma -\frac{x^2}{2})\end{array}$$

No wait, that’s wrong, we have $\tau$ on both sides of the update equation. So in fact, this ain’t conjugate. I wonder what is?