R/forwardAlg.R

In jiangrongo/HMM:

#' forwardAlg
#'
#' \code{forwardAlg} performs the Forward Algorithm on a single trajecture of observed chain X
#'
#' \code{forwardAlg} computes the matrix, alpha, where alpha[k,t] = P(Z[t]=k, X[1:t]=x[1:t])
#'
#' @export
#' @param X     a vector of observed states
#' @param trans a matrix of transition probability
#' @param u     a vector of means of Normal distribution for each state (emission probability)
#' @param sig   a vector of standard deviations of Normal distribution for each state (emission probability)
#' @return      A matrix, alpha.
#'
#' @examples
#' set.seed(1221)
#' df <- generateHMM(num=2,n=10)
#' forwardAlg(df$X[,1], df$trans, df$u, df$sig)

forwardAlg <- function(X, trans, u, sig){
  n <- length(X)     # length of Markov chain
  M <- dim(trans)[1] # number of states of Markov chain

  alpha <- matrix(0, M, n)
  alpha[,1] <- dnorm(X[1], u, sig)/M

  # compute alpha by time/column
  for (t in 2:n){
    Pxz <- dnorm(X[t], u, sig)
    Pzza <- t(trans) %*% alpha[,t-1]
    alpha[,t] <- Pxz*Pzza
  }

  return(alpha)
}