R/dQ_em.R

In ick003/SpTMixture: Mixture models with spatio-temporal dependencies

#' Derivative function of the Q function of an EM algorithm
#'
#' Derivative function of the Q function for a mixture (EM)
#' @param theta parameters of the mixture that need estimation
#' @param y observations
#' @param p.theta result from the expectation step (EM - probabilities)
#' @param mean_func function to calculate the mean for each mixture. The result is a matrix where the number of columns is the number of element in the mixture.
#' @keywords mixture model
#' @references Dempster, A. and Laird, N. and Rubin, D. (1977) Maximum likelihood from incomplete data via the EM algorithm. \emph{JRSS B} \bold{39} 1--38.
#' @export
#' @examples
#' dQ_em()
dQ_em = function(theta,y, p.theta, mean_func, ...){

  ID = levels(y$ID)

  m = mean_func(theta, ...)

  sigma = theta[(length(theta) - ncol(p.theta) + 1):length(theta)]

  dq = rep(0, length(theta))
  for(i in 1:nrow(p.theta)){
    y.p = y[y$ID == ID[i], ]
    for(j in 1:ncol(p.theta)){

      idx.j = (1:(length(theta) / ncol(p.theta) - ncol(p.theta))-1)*(length(theta) / ncol(p.theta) - ncol(p.theta)) +  j
      m.p = m[, j]
      Xj =
      dobs = -2 * t(Xj) %*% (y.p$obs - m.p) + t(y.p$obs - m.p) %*% Xj

      dq[idx.j] = dq[idx.j] + p.theta[i,j] / sigma[j] * dobs

    }

  }
  if(is.na(q)){q = -Inf}
  return(q)
}