HWPerm: Permutation test for Hardy-Weinberg equilibrium

In HardyWeinberg: Statistical Tests and Graphics for Hardy-Weinberg Equilibrium

Permutation test for Hardy-Weinberg equilibrium

Description

Function HWPerm does a permutation test for Hardy-Weinberg equilibrium using a user-supplied test statistic.

Usage

HWPerm(x, nperm = 17000, verbose = TRUE, x.linked = FALSE,
FUN = ifelse(x.linked,Chisquare.x,Chisquare), eps=1e-10, ...) 

Arguments

x	A vector of genotype counts (AA,AB,BB)
nperm	The number of permutations
verbose	verbose = TRUE will print results, verbose = FALSE is silent.
x.linked	x.linked=FALSE indicates the marker is autosomal (default), and x.linked=TRUE indicates it resides on the X-chromosome.
FUN	An function call for calculating the test statistic for HWE (see examples below)
eps	Tolerance for comparison of floating point numbers (1e-10 by default)
...	Additional parameters for the function call argument FUN

Details

The set of alleles for the observed sample is permuted. Consequently, the test is conditional on allele frequency.

Value

HWPerm returns a list with the components:

stat	value of the chosen test statistic for the observed sample.
pval	p-value of the permutation test.

Author(s)

Jan Graffelman [e-mail]

References

Ziegler, A. & Konig, I.R. (2006) A statistical approach to genetic epidemiology. Wiley.

See Also

HWChisq,HWExact,HWLratio

Examples

x <- c(MM=298,MN=489,NN=213)
## Not run: 
#
# use a default chi-square statistic
#
HW.test <- HWPerm(x,nperm=10000,verbose=TRUE)
#
# use a chi-square statistic with continuity correction.
#
HW.test <- HWPerm(x,nperm=10000,verbose=TRUE,
FUN=function(z) HWChisq(z,verbose=FALSE)$chisq,cc=0.5) # 
#
# use a likelihood ratio statistic.
#
HW.test <- HWPerm(x,nperm=10000,verbose=TRUE,
FUN=function(y) HWLratio(y,verbose=FALSE)$G2)
# 
# use an exact test p-value
#
HWPerm(x,nperm=10000,verbose=TRUE,FUN=function(y) 1-HWExact(y,verbose=FALSE)$pval)
# 
#
# Permutation test for a marker on the X chromosome
#
rs5968922 <-  c(A=392, B=212, AA=275, AB=296, BB=80)
HW.test <- HWPerm(rs5968922,nperm=10000,x.linked=TRUE,verbose=TRUE)

## End(Not run)

HardyWeinberg documentation built on May 29, 2024, 6:17 a.m.