gibbs_gl {hwep}

Gibbs sampler under random mating using genotype log-likelihoods.

Description

Gibbs sampler under random mating using genotype log-likelihoods.

Usage

gibbs_gl(
  gl,
  alpha,
  B = 10000L,
  T = 1000L,
  more = FALSE,
  lg = FALSE,
  verbose = TRUE
)

Arguments

gl

The matrix of genotype log-likelihoods. The columns index the dosages and the rows index the individuals. gl[i,j] is the genotype log-likelihood for individual i at dosage j. It is assumed that natural log is used.

alpha

Vector of hyperparameters for the gamete frequencies. Should be length (x.length() - 1) / 2 + 1.

B

The number of sampling iterations.

T

The number of burn-in iterations.

more

A logical. Should we also return posterior draws (TRUE) or not (FALSE).

lg

Should we return the log marginal likelihood (true) or not (false).

verbose

A logical. Should we print the progress?

Value

A list with some or all of the following elements

• mx: The estimate of the marginal likelihood

• p_tilde: The value of p used to evaluate the posterior density.

• p: The samples of the gamete frequencies

• z: The samples of the individual genotypes

• post: The samples of the full conditionals of p_tilde.

Author(s)

David Gerard

Examples

set.seed(1)
ploidy <- 8

## Simulate under the null
p <- stats::runif(ploidy / 2 + 1)
p <- p / sum(p)
q <- stats::convolve(p, rev(p), type = "open")
nvec <- c(stats::rmultinom(n = 1, size = 100, prob = q))
gl <- simgl(nvec)

gibbs_gl(gl = gl, alpha = rep(1, ploidy / 2 + 1), lg = TRUE)