How to calculate log likelihood for gaussian mixture model I'm trying to check the calculation of the log likelihood of a 2 component Gaussian Mixture Model using optim, but I get the wrong answer (it should return mu, sigma, alpha actual). 
The log likelihood is calculated as: 
$log(p(y_i|\alpha,\mu,\sigma))=\sum_{i=1}^Nlog(\sum_{k=1}^2\alpha_k*N(y_i|\mu_k,\sigma_k))$
using the method shown on this blog. Thanks for any help.
#generate response
mu_actual <- c(-5,3)
sigma_actual <- c(2,3)
alpha_actual <- c(0.3,0.7)

comp <- rbinom(n=1e3,size=1,prob=alpha_actual[1]) 
Y <- rnorm(1e3, mean=mu_actual[comp+1], sd=sigma_actual[comp+1])

#likelihood calculation
fn_likelihood_optim <- function(theta,Y){
  mu_1 <- theta[1]
  mu_2 <- theta[2]
  sigma_1 <- theta[3]
  sigma_2 <- theta[4]
  alpha_1 <- theta[5]
  alpha_2 <- 1-theta[5]

  like_1 <- log(alpha_1)+dnorm(Y,mu_1,sigma_1,log=TRUE)
  like_2 <- log(alpha_2)+dnorm(Y,mu_2,sigma_2,log=TRUE)
  log_like <- sum(log(exp(like_1)+exp(like_2)))

  if(is.infinite(log_like) & log_like>0) log_like <- 1e6
  if(is.infinite(log_like) & log_like<0) log_like <- 1e-6

  return(-log_like)
}

#check likelihood with optim - should return theta actual
z <- optim(par = c(0,5,1,2,0.5),
           fn = fn_likelihood_optim,
           method = "L-BFGS-B",
           lower=c(-10,-10,0.001,0.001,0),
           upper=c(10,10,10,10,0.99),
           Y=Y)

z$par

 A: You are getting infinity if alpha_1 every touches zero, so set the alpha_1 lower boundary to 0.01 or something not zero. I think there is a typo here if(is.infinite(log_like) & log_like>0) log_like <- 1e6, should be 1e-6.
set.seed(111)
comp <- rbinom(n=1e3,size=1,prob=alpha_actual[1]) 
Y <- rnorm(1e3, mean=mu_actual[comp+1], sd=sigma_actual[comp+1])

fn_likelihood_optim <- function(theta,Y){
  mu_1 <- theta[1]
  mu_2 <- theta[2]
  sigma_1 <- theta[3]
  sigma_2 <- theta[4]
  alpha_1 <- theta[5]
  alpha_2 <- 1-theta[5]

  like_1 <- log(alpha_1)+dnorm(Y,mu_1,sigma_1,log=TRUE)
  like_2 <- log(alpha_2)+dnorm(Y,mu_2,sigma_2,log=TRUE)
  log_like <- sum(log(exp(like_1)+exp(like_2)))

  if(is.infinite(log_like) & log_like>0) log_like <- 1e-6
  if(is.infinite(log_like) & log_like<0) log_like <- 1e-6

  return(-log_like)
}

#check likelihood with optim - should return theta actual
z <- optim(par = c(0,5,1,2,0.5),
           fn = fn_likelihood_optim,
           method = "L-BFGS-B",
           lower=c(-10,-10,0.001,0.001,0.01),
           upper=c(10,10,10,10,0.99),
           Y=Y)

z$par
[1]  3.0831575 -5.0493771  2.7413426  1.9876056  0.2740447

