EM Algorithm Negative Log Likelihood

I am fitting a shared parameter model where I have a binary $$Y_2$$ and a continuous $$Y_1$$, and I am using random effects in the mean structure to model them jointly

$$f(Y_1,Y_2)=\int f(Y_1,Y_2,b)\,\mathrm d b=\int f(Y_1|b)f(Y_2|b)f(b)\,\mathrm d b$$

where the mean structure for each

$$E[Y_1|b]=X^t\beta +b$$ $$E[Y_2|b]=\Phi(X^t\beta+b)$$

where $$Y_1$$ is normal and $$Y_2$$ is Bernoulli.

I am using the EM algorithm to solve for this, and it is increasing the $$\log$$ likelihood monotonically but the log-likelihood is very negative $$(-20100)$$ in the beginning and is increasing slowly. I am using the marginal linear and glm coefficient estimates (not considering the shared random effects b) as starting points.

What could be causing the large negative likelihood and if it converges at a negative log likelihood is there something wrong with my model?

• How are you modeling $f(b)$? – A rural reader Mar 1 at 1:42
• It is being modeled as a normal random variable. – Noe Vidales Mar 1 at 1:43
• Thanks! Maybe one more thing: what does your Bernoulli $f(y_2 \, |\, b)$ look like? I see your model for the mean, but the standard Bernoulli has mean between zero and one; at the moment these seem inconsistent, any light will help. 🙂 – A rural reader Mar 1 at 1:55
• Sorry forgot to include probit link – Noe Vidales Mar 1 at 2:13
• One more thing. Is it correct to infer that for each sample you have measurements of $x_1, \dots, x_m$ and $y_1$ and $y_2$? Then for $y_1$, for example, you're modeling it as normally distributed with mean $\beta_1 x_1 + \cdots \beta_m x_m + b$ and variance $\sigma^2$. Same sort of thing for the Bernoulli. So the parameters you're estimating include the $\beta$s, $\sigma^2$. – A rural reader Mar 1 at 2:44

A complete log-likelihood would look like \begin{align*}\ell^*(\beta) &= -\frac{1}{2\sigma^2}\sum_{i=1}^n (y_{1i}-x_i^\text{T}\beta-b_i)^2\\ &\ -\frac{1}{2}\sum_{i=1}^n (y^*_{2i}-x_i^\text{T}\beta-b_i)^2\\ &\ -\frac{1}{2\tau^2}\sum_{i=1}^n b_i^2-\frac{n}{2}\log \tau^2-\frac{n}{2}\log \sigma^2\end{align*} with $$(i=1,\ldots,n)$$ $$y_{2i}=\mathbb I_{y_{2i}^*\ge 0}$$ and hence the E-step involves \begin{align*}Q(\beta,\sigma,\tau)&=n\log \sigma^2+n\log \tau^2+ \mathbb E\left[\frac{1}{\tau^2}\sum_{i=1}^n B_i^2 \Big|\mathbf y_1,\mathbf y_2\right]\\ &+\mathbb E\left[\frac{1}{\sigma^2}\sum_{i=1}^n (y_{1i}-x_i^\text{T}\beta-B_i)^2\Big|\mathbf y_1,\mathbf y_2\right]\\ &+\mathbb E\left[\sum_{i=1}^n (Y^*_{2i}-x_i^\text{T}\beta-B_i)^2\Big|\mathbf y_1,\mathbf y_2\right]\\ \end{align*} Note that there is not absolute scale or range for a log-likelihood, hence its numerical value being "large" or negative does not mean anything per se. In fact, a normal -log likelihood such as $$-\sum_{i=1}^n (X_i-\mu)^2/\sigma^2-n\log(2\pi\sigma^2) \qquad X_i\sim\mathcal N(\mu,\sigma^2)$$ is of order $$O(-n)$$ hence goes almost surely to $$-\infty$$ as $$n$$ increases.
To check whether or not the EM algorithm is acting correctly, you should work with a simulated sample from the model (with known parameters $$\beta,\sigma,\tau$$) so that you can compute the complete log-likelihood value at the true parameter value (which gives the scale of the likelihood at the maximum) and observe if EM converges to nearby values.