Past exam of the mathematics course of the University of Cambridge / 2002 / iii / Paper 40 / 4

Put $v={\sigma }^2>0$, $\bar{x}=n^{-1}{\sum }_i{x}_i$ and $S={\sum }_i(x_i-\bar{x}{)}^2$. The likelihood function and log-likelihood areFor every fixed $v$, the unique maximum in $\mu$ is $\bar{x}$. At that value, $\partial \ell /\partial v=-n/(2v)+S/(2v^2)$ changes sign from positive to negative at $v=S/n$. Thus, for $S>0$,This is the global maximum-likelihood estimate, since the first maximization was global for each $v$ and the second has a unique peak. If $S=0$, the likelihood function is unbounded at $\mu =\bar{x}$, $v\downarrow 0$; there is no maximizer in the parameter space $v>0$. The freezing argument below uses the nondegenerate case $S>0$.

At temperature $T>0$, target a probability density function proportional to $\exp (f(\theta )/T)$ on a parameter region where it is normalizable, possibly with a specified proper base measure. A Metropolis–Hastings algorithm proposal $q(\theta ,{\theta }^{\prime })$ has acceptance probabilityFor a symmetric proposal, improving moves are accepted and a loss $\Delta$ is accepted with probability $e^{-\Delta /T}$. Run updates while gradually lowering $T$ toward zero. This simulated annealing scheme initially permits exploration and increasingly favors high objective values. Appropriate cooling and mixing conditions are needed for a general global convergence theorem; a finite arbitrary schedule does not guarantee finding a global maximum.

To maximize a positive likelihood function, take $f=\log L$, which has the same maximizers. The resulting likelihood-power annealing target is $L^{1/T}$, and available Gibbs sampler conditional distributions can replace Metropolis proposals.

Use the Gaussian likelihood Gibbs annealing target ${\pi }_T(\mu ,v)\propto L(\mu ,v{)}^{1/T}$ relative to $d\mu dv$, with $S>0$. Completing the normal square and recognizing an inverse-gamma kernel giveThe shape is $n/(2T)-1$, because the probability density function in $v$ is $v^{-n/(2T)}{e}^{-[S+n(\mu -\bar{x}{)}^2]/(2Tv)}$. The base measure matters: inserting a prior factor $1/v$ would change this shape.

Integrating out $\mu$ contributes a factor proportional to $v^{1/2}$, so the marginal is $v\sim \mathrm{IG}(n/(2T)-3/2,S/(2T))$. Thus the joint target is proper for $T<n/3$. Choose a decreasing sequence $0<T_m\le n/10$ tending to zero and a finite starting value $v_0>0$. At iteration $m$, sample ${\mu }_m$ from the first conditional using $v_{m-1}$, then sample $v_m$ from the second using ${\mu }_m$. This is an explicit annealing algorithm entirely using Gibbs sampling, with all its inverse-gamma moments finite for the convergence proof.

We prove convergence of the actual sequential updates above, not just concentration of ideal equilibrium distributions. Write $M_m=\mathbb{E}{v}_m$ and $N_m=\mathbb{E}{v}_m^2$. Given $v_{m-1}$, the centered normal draw satisfiesThe conditional inverse-gamma distribution moments therefore giveFor $T_m\le n/10$, the coefficient multiplying $M_{m-1}$ is at most $1/6$, so $M_m$ stays bounded by a geometric-recursion argument. The coefficient multiplying $N_{m-1}$ is at most $1/8$, and the other terms are bounded using the first-moment bound, so $N_m$ also stays bounded. These bounds start from the finite deterministic $v_0$.

Since $T_m\to 0$, the boundedness makes every term containing $T_m{M}_{m-1}$ or $T_m^2{N}_{m-1}$ vanish. Hence $M_m\to S/n$ and $N_m\to S^2/n^2$. ConsequentlyThusThis establishes mean-square convergence for the specified Gaussian annealer. It does not assert that every implementation or cooling rule for every objective has the same guarantee, nor that mean-square convergence alone is almost-sure convergence.

For observed data $x$ and latent variables $Z$, let ${\ell }_c(\theta ;x,Z)=\log p_{\theta }(x,Z)$ be the complete-data log-likelihood. From a current parameter ${\theta }^{(m)}$, the expectation-maximization algorithm computesand then chooses ${\theta }^{(m+1)}$ maximizing this function in its first argument. Thus the E-step averages the missing information conditional on the current model and observed data; the M-step treats those expectations as fixed while maximizing.

For a short proof of likelihood function monotonicity, let $q_m(Z)=p_{{\theta }^{(m)}}(Z\mid x)$. ThenThe Jensen inequality for the logarithm gives $\ell (\theta ;x)-\ell ({\theta }^{(m)};x)\ge Q(\theta \mid {\theta }^{(m)})-Q({\theta }^{(m)}\mid {\theta }^{(m)})$. The M-step makes the right side nonnegative. Thus observed likelihood function cannot decrease. This proves the usual EM likelihood monotonicity, while allowing convergence to local stationary points rather than guaranteeing a global maximum.

Put ${\alpha }_1=\alpha$, ${\alpha }_2=1-\alpha$, and introduce a categorical class label for each observation, encoded by $Z_{i1},Z_{i2}\in \{0,1\}$ with $Z_{i1}+Z_{i2}=1$. Its prior class probability is ${\alpha }_j$, and conditional on class $j$ the observation probability density function is $f_j(x_i)$. Thus its joint observation-label probability density function isIndependence across observations yields the finite Gaussian mixture with a common variance complete-data likelihood functionSumming over the two possible labels at each observation recovers ${\prod }_i[\alpha f_1(x_i)+(1-\alpha )f_2(x_i)]$, the observed-data likelihood function. The product requested in the PDF is therefore a joint complete-data likelihood function. Strictly conditioning on $Z$ would remove the factors ${\alpha }_j$; the PDF's conditional notation for this joint product is imprecise.

The E-step uses the posterior class probabilitiesLet $n_j^{(m)}={\sum }_i{w}_{ij}^{(m)}$, so $n_1^{(m)}+n_2^{(m)}=n$. Up to constants, the expected complete-data log-likelihood isMaximizing the first term with ${\alpha }_1+{\alpha }_2=1$ gives ${\alpha }^{(m+1)}=n_1^{(m)}/n$. Differentiating in each ${\mu }_j$ gives ${\mu }_j^{(m+1)}={\sum }_i{w}_{ij}^{(m)}{x}_i/n_j^{(m)}$. After these updates, maximizing in the shared variance givesThe residuals use the new component means, while the weights remain from the old parameters. The divisor is $n$, not $n-2$, because this is a maximum-likelihood update. Starting with $0<\alpha <1$ and $v>0$ gives positive responsibilities and positive $n_j$; a zero effective component size at a boundary would leave its mean unidentified. These are the EM for Gaussian mixtures with a common variance transitions.