R/forward_prob.R

Defines functions cat_Mult_HMM_fw

#' @keywords internal
# Could maybe made external
# Calculates the forward probabilities, used for sampling the state sequence
# Based on Zuchini 2016.

cat_Mult_HMM_fw <- function(x, m, emiss, n_dep, gamma, delta = NULL){
  if(is.null(delta)) {
    delta <- solve(t(diag(m) - gamma + 1), rep(1, m))
  }
  n        <- dim(x)[1]
  lalpha   <- alpha_prob <- matrix(NA, m, n)
  allprobs <- matrix(NA, nrow = n, ncol = m)
  inp      <- matrix(ncol = n_dep, nrow = m)
  for(i in 1:n){
    for (q in 1:n_dep){
      inp[, q]    <- emiss[[q]][, x[i, q]]
    }
    allprobs[i, ] <- apply(inp, 1, prod)
  }
  foo             <- delta * allprobs[1, ]
  sumfoo          <- sum(foo)
  alpha_prob[, 1] <- foo/sumfoo
  lscale          <- log(sumfoo)
  lalpha[, 1]     <- log(alpha_prob[, 1]) + lscale
  for (i in 2:n){
    foo              <- alpha_prob[, (i - 1)] %*% gamma * allprobs[i, ]
    sumfoo           <- sum(foo)
    alpha_prob[, i]  <- foo / sumfoo
    lscale           <- lscale + log(sumfoo)
    lalpha[, i]       <- log(alpha_prob[, i]) + lscale
  }
  list(la = lalpha, forward_p = alpha_prob)
}

Try the mHMMbayes package in your browser

Any scripts or data that you put into this service are public.

mHMMbayes documentation built on Oct. 30, 2019, 5:05 p.m.