Simulate samples from a Dirichlet prior or posterior under HWE

Share:

Description

Function to simulate samples from the HWE Dirichlet model. Can be used for samples from the prior or the (conjugate) Dirichlet posterior, both in the k allele case. Samples are generated for the allele frequencies in the order p_{1},p_{2},...,p_{k}.

Usage

1
DirichSampHWE(nvec, bvec0, nsim)

Arguments

nvec

vector of genotype frequencies in the order n_{11}, n_{12},..., n_{1k},n_{22} ..., n_{2k},..., n_{kk}.

bvec0

vector of length k Dirichlet prior parameters, where k is the number of alleles.

nsim

number of samples to simulate from the prior/posterior.

Details

Uses the rdirichlet function from the MCMCpack library.

Value

pvec

matrix of size nsim \times k containing samples for the genotype frequencies, in the order p_{1}, p_{12},..., p_{k}.

Author(s)

Jon Wakefield (jonno@u.washington).

References

Wakefield, J. (2010). Bayesian methods for examining Hardy-Weinberg equilibrium. Biometrics; Vol 66:257-65

See Also

DirichSampSat, DirichNormSat, DirichNormHWE

Examples

1
2
3
4
5
6
7
8
9
# First sample from the prior
PriorSampHWE <- DirichSampHWE(nvec=rep(0,10),bvec0=rep(1,4),nsim=1000)
par(mfrow=c(1,1))
hist(PriorSampHWE$pvec[,1],xlab="p1",main="")
# Now sample from the posterior
data(DiabRecess)
PostSampHWE <- DirichSampHWE(nvec=DiabRecess,bvec0=rep(1,4),nsim=1000)
par(mfrow=c(1,1))
hist(PostSampHWE$pvec[,1],xlab="p1",main="")