forward: Applies the forward algorithm to a genotype dataset

In EpiSlim/epiGWAS: Robust Methods for Epistasis Detection

Description

Applies the forward_sample function to each row in X. If the ncores > 1, the function calling is performed in a parallel fashion to reduce the running time. The parallelization backend is doParallel. If the latter package is not installed, the function switches back to single-core mode.

Usage

forward(X, p_init, p_trans, p_emit, ncores = 1)

Arguments

X	genotype matrix. Each row corresponds to a separate sample
p_init	marginal distributions for the first hidden state
p_trans	3D dimensional array for the transition probabilities
p_emit	3D dimensional array for the emission probabilities
ncores	number of threads (default 1)

Value

A vector of log probabilities

References

Rabiner, Lawrence R. 'A tutorial on hidden Markov models and selected applications in speech recognition.' Proceedings of the IEEE 77.2 (1989): 257-286.

Examples

p <- 3 # Number of states
K <- 2 # Dimensionality of the latent space

p_init <- rep(1 / K, K)
p_trans <- array(runif((p - 1) * K * K), c(p - 1, K, K))
# Normalizing the transition probabilities
for (j in seq_len(p - 1)) {
  p_trans[j, , ] <- p_trans[j, , ] / (matrix(rowSums(p_trans[j, , ]), ncol = 1) %*% rep(1, K))
}

p_emit <- array(stats::runif(p * 3 * K), c(p, 3, K))
# Normalizing the emission probabilities
for (j in seq_len(p)) {
  p_emit[j, , ] <- p_emit[j, , ] / (matrix(rep(1, 3), ncol = 1) %*% colSums(p_emit[j, , ]))
}

n <- 2
X <- matrix((runif(n * p, min = 0, max = 1) < 0.4) +
            (runif(n * p, min = 0, max = 1) < 0.4), nrow = 2)

# Computing the joint log-probabilities
log_prob <- forward(X, p_init, p_trans, p_emit)

EpiSlim/epiGWAS documentation built on Nov. 19, 2019, 7:15 p.m.