HWEImportSamp: Importance sampling to calculate the normalizing constant...

In HWEBayes: Bayesian investigation of Hardy-Weinberg Equilibrium via estimation and testing.

Description

Importance sampling to calculate the normalizing constant under the single f model. Two proposals are available, either sampling from the prior or sampling from a normal distribution whose mean vector and variance-covariance matrix must be specified. The latter may be taken from an MCMC analysis using, for example, WinBUGS. In all cases the likelihood is multinomial and the prior is Dirichlet on the allele frequencies, and normal on λ where λ=\log(f-f_{\min})/(1-f)). See Weir (1996) for a description of HWE and different models/parameterizations.

Usage

HWEImportSamp(nsim, nvec, ischoice, lambdamu, lambdasd, alpha,
gmu = rep(0, length(alpha)), gsigma = diag(0, nrow = length(alpha), 
ncol = length(alpha)))

Arguments

nsim	the number of points to sample to calculate the estimate.
nvec	vector of genotype frequencies in the order n_{11}, n_{12},..., n_{1k},n_{22} ..., n_{2k},..., n_{kk}.
ischoice	choice of importance sampling proposal, =1 gives a normal distribution with mean and variance that must be specified (as gmu and gsigma) and =2 is from the prior.
lambdamu	the mean of the prior for λ.
lambdasd	the variance of the prior for λ.
alpha	the vector of k parameters for the Dirichlet prior on the allele frequencies.
gmu	the mean of the importance sampling proposal, of length k, where k is the number of alleles.
gsigma	the variance of the importance sampling proposal, a matrix of dimension k \times k, where k is the number of alleles.

Value

PrnH1	the estimate of the normalizing constant
varest	the variance of the estimate of the normalizing constant

Warning

As always with importance sampling the procedure can be very unstable, particularly for large k. Hence rerunning the function with different simulation sample sizes, and different gmu and gsigma is recommended

Author(s)

Jon Wakefield ([e-mail])

References

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

Weir, B.S. (1996). Genetic Data Analysis II. Sunderland MA: Sinauer.

See Also

LambdaOptim, DirichNormSat, DirichNormHWE, TriangNormHWE

Examples

alpha <- c(1,1,1,1) # prior on allele frequencies
# gmu and gsigma were obtained from a WinBUGS run
gmu <- c(-0.4633092,0.3391625,0.3397936,-3.5438008)
gsigma <- matrix(c(
0.07937341,0.02819656,0.02766583,0.04607996,
0.02819656,0.07091320,0.04023827,0.01657028,
0.02766583,0.04023827,0.07042278,0.01752266,
0.04607996,0.01657028,0.01752266,0.57273683),nrow=4,ncol=4)
data(DiabRecess)
HWEImportSamp(nsim=5000,nvec=DiabRecess,ischoice=1,lambdamu=-2.95,
   lambdasd=1.07,alpha=alpha,gmu,gsigma)
HWEImportSamp(nsim=5000,nvec=DiabRecess,ischoice=2,lambdamu=-2.95,
   lambdasd=1.07,alpha=alpha)

HWEBayes documentation built on May 2, 2019, 11:07 a.m.