cond_prob: Computes the propensity scores
In epiGWAS: Robust Methods for Epistasis Detection

Description Usage Arguments Value See Also Examples

In this function, and for each sample, we compute both propensity scores P(A=1|X) and P(A=0|X). The application of the forward algorithm on the passed hmm allows us to estimate the joint probability of (A, X), for all values of the target variant A = 0, 1, 2. The Bayes formula yields the corresponding conditional probabilities. Depending on the binarization rule, we combine them to obtain the propensity scores.

1	cond_prob(X, target_name, hmm, binary = FALSE, ncores = 1)

`X`	genotype matrix. Make sure to assign `colnames(X)` beforehand.
`target_name`	target variant name
`hmm`	fitted parameters of the fastPHASE hidden Markov model. The HMM model is to be fitted with the `fast_HMM` function.
`binary`	if `TRUE`, the target SNP values 0 and (1,2) are respectively mapped to 0 and 1. That describes a dominant mechanism. Otherwise, if `FALSE`, we encode a recessive mechanism where the values 0 and 1 respectively map to (0,1) and 2.
`ncores`	number of threads (default 1)

Two-column propensity score matrix. The first column lists the propensity score P(A=0|X), while the second gives P(A=1|X).

fast_HMM

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, , ]))
}

hmm <- list(pInit = p_init, Q = p_trans, pEmit = p_emit)

n <- 2
X <- matrix((runif(n * p, min = 0, max = 1) < 0.4) +
            (runif(n * p, min = 0, max = 1) < 0.4),
            nrow = 2, dimnames = list(NULL, paste0("SNP_", seq_len(p))))

cond_prob(X, "SNP_2", hmm, ncores = 1, binary = TRUE)