# Calculate log likelihood of mixture of gaussians “by hand” in R

I'm trying to ensure the the method of calculating log likelihood for a model produced using mixtools vs a model produced using MLE estimates of mu and sigma are the same. The best way I can think of doing this is to write my own function to calculate log likelihood.

The first model, produced by mixtools' normalmixEM function, is set to be a mixture of two gaussians. I'm struggling to find how log likelihood is calculated in the case of mixtures, taking into account the lambda values of each gaussian (i.e. its "contribution").

The second model is set to be a single gaussian whose parameters are estimated by using simple MLE estimates for mu and sigma based on the data.

I've seen it suggested that the likelihood can be calculated for a single gaussian using the call:

ll=sum(dnorm(y,mu,sigma,log=T))


I'm not sure how this call is the log likelihood, given I know about how LL is calculated.

Can someone confirm whether this actually gives me the LL, and if so, how I can integrate this method of calculating LL with a mixture of two gaussians?

For the log-likelihood of a mixture of Gaussians, you may want to have a look here:

https://stephens999.github.io/fiveMinuteStats/intro_to_em.html

sum(dnorm(y, mu, sigma, log=T)) actually returns the log-likelihood.

For a mixture of 2 gaussian you can use the function

compute.log.lik <- function(L, w) {
L[,1] = L[,1]*w[1]
L[,2] = L[,2]*w[2]
return(sum(log(rowSums(L))))
}


where for instance

L = matrix(NA, nrow=length(X), ncol= 2)
L[, 1] = dnorm(X, mean=mu.true[1], sd = sigma.true[1])
L[, 2] = dnorm(X, mean=mu.true[2], sd = sigma.true[2])

• Thank you! Can you perhaps elaborate on why 'sum(dnorm(y, mu, sigma, log=T))' is the LL? Is each new data point in dnorm calculated by taking the difference between each data point and mu? – Rebecca Eliscu Jun 23 '20 at 19:03
• It's basically the sum of the log norm that comes from the log of the product of iid normally distributed points. par(mfrow = c(2,1)) plot(function(x) dnorm(x, log = TRUE), -60, 50, main = "log { Normal density }") curve(log(dnorm(x)), add = TRUE, col = "red", lwd = 2) shows the equivalence between log(dnorm()) and dnorm(..., log=T) – user289381 Jun 23 '20 at 19:11
• I see. Apologies, I realize now that my confusion stems from a fundamental lapse in my memory about how a Gaussian is defined to begin with (with each data point subtracted by mu in the exp() part of the function). Thanks so much for the help! – Rebecca Eliscu Jun 23 '20 at 19:23
• May I ask you to mark my answer as useful, please? :) – user289381 Jun 23 '20 at 19:52
• I did, but unfortunately: "Thanks for the feedback! Votes cast by those with less than 15 reputation are recorded, but do not change the publicly displayed post score." – Rebecca Eliscu Jun 23 '20 at 20:19