LD.clump {gaston}

LD clumping

Description

Construct group of SNPs in LD with 'top associated SNPs'

Usage

LD.clump(x, p, r2.threshold, p.threshold, max.dist = 500e3)

Arguments

x

A bed.matrix

p

A vector of p-values, or a data frame including p-values, such as sent back by association.test

r2.threshold

The maximum LD (measured by r^2) between SNPs in a group

p.threshold

The threshold used to define associated SNPs

max.dist

The maximum distance for which the LD is computed

Details

The p-values provided through argument p are assumed to correspond to the result of an association test with the SNPs of x.

The aim of the function is to construct cluster of SNPs in strong LD with associated SNPs. The algorithm first seeks the SNP with the lowest p-value (below p.threshold) ; this SNP will be the 'index' of a cluster. The corresponding cluster is constructed by aggregating SNPs that are in LD (above r2.threshold) with the index. The cluster's name is the position of the index SNP. The processus is repeated on the SNPs which are not yet attributed to a cluster, until there is no associated SNP (ie SNP with a p-value below threshold) left. The remaining SNPs are attributed to cluster 0.

The LD is computed only for SNP pairs for which distance is inferior to max.dist, expressed in number of bases: above this distance it is assumed to be null.

Value

If p was a data frame, then the function returns the same data frame with to extra columns, cluster and is.index. If p was a vector of p-values, it returns a data frame with columns chr, id, pos, p, cluster and is.index.

See Also

LD, LD.thin

Examples

# Construct a bed matrix
x <- as.bed.matrix(TTN.gen, TTN.fam, TTN.bim)
standardize(x) <- "p"
     
# simulate quantitative phenotype with effect of SNPs #108 and #631
beta <- numeric(ncol(x))
beta[c(108,631)] <- 0.5
set.seed(1)
y <- x %*% beta + rnorm(nrow(x))
     
# association test with linear model 
test <- association.test(x, y, method="lm", response = "quanti")

# LD clumping
test <- LD.clump(x, test, r2.threshold = 0.25, p.threshold = 1e-8)

# use as.factor for a quick-and-dirty cluster colouring on the manhattan plot 
manhattan(test, col = as.factor(test$cluster), pch = 20)