Skip to contents

GBS/RADseq short and long distance linkage disequilibrium pruning

Source: R/filter_ld.R
filter_ld.Rd

SNP short and long distance linkage disequilibrium pruning.

Filter target: Markers.

What sets appart radr LD pruning is the RADseq data tailored arguments:

  • minimize short linkage disequilibrium (LD): 5 values available for filter.short.ld argument (see below).

  • reduce long distance LD: Long distance LD pruning is usually advised to avoid capturing the variance LD in PCA analysis.

    Use the argument filter.long.ld with values between 0.7 and 0.9 is a good starting point. Ideally, you want to visualize the LD before choosing a threshold.

    Strategically, run the function with filter.long.ld argument and refilter the data using the outlier statistic generated by the function (printed on the figure in the output) and using long.ld.missing = TRUE. This advanced argument will choose the best SNP based on missing data statistics, instead of choosing randomly one SNP (see details).

This function is used internally in radr and might be of interest for users.

Usage

filter_ld(
  data,
  interactive.filter = TRUE,
  filter.short.ld = "mac",
  filter.long.ld = NULL,
  parallel.core = parallel::detectCores() - 1,
  filename = NULL,
  verbose = TRUE,
  ...
)

Arguments

data

A tidy genomic data frame or another genomic object supported by the calling function.

interactive.filter

Logical indicating whether an interactive filtering session may display diagnostics and ask for thresholds. Default: interactive.filter = TRUE.

filter.short.ld

(character) Five options:

  1. filter.short.ld = "random" for a random selection of 1 SNP on the read,

  2. filter.short.ld = "first" for the first one on the read...,

  3. filter.short.ld = "last" for the last SNP on the read and

  4. filter.short.ld = "middle" for locus with > 2 SNPs/read the option to select at random one SNP between the first and the last SNP on the read. If the locus as <= 2 SNPs on the read, the first one is selected. Note that for that last option, the numbers are reported.

  5. filter.short.ld = "mac" will select the SNP on the locus with the maximum global Minor Allele Count (MAC).

    Using filter.short.ld = NULL, skip this filter.

Default: filter.short.ld = "mac".

filter.long.ld

(optional, double) The threshold to prune SNP based on Long Distance Linkage Disequilibrium. The argument filter.long.ld is the absolute value of measurement. Default: filter.long.ld = NULL.

parallel.core

Number of workers available for parallel operations. Default: parallel.core = parallel::detectCores() - 1.

filename

(optional, character) File name prefix for file written in the working directory. Default: filename = NULL.

verbose

Logical indicating whether progress messages are emitted. Default: verbose = TRUE.

...

Additional arguments passed to lower-level screening or filtering functions.

Value

The filtered data in the same representation as the input. GDS marker metadata and active variants are updated in place. LD summaries, plots, marker lists, and filtering parameters are written to the output folder when applicable.

Details

The function requires SNPRelate (see example below on how to install).

Advanced mode, using dots-dots-dots

  • maf.data (path) this argument is no longer supported. It's a small cost in time in favour of making sure the MAC/MAF fits the actual data.

  • long.ld.missing (logical) With long.ld.missing = TRUE. The function first generates long distance LD values between markers along the same chromosome or scaffold with SNPRelate::snpgdsLDMat. Based on the LD threshold (filter.long.ld) SNPs in LD will be pruned based on missingness. e.g. if 4 SNPs are in LD, the 1 SNP selected in the end is base on genotyping rate/missingness. If this statistic is equal between the SNPs in LD, 1 SNP is chosen randomly.

    Using missigness add extra computational time. To speed the analysis when missingness between markers is not an issue, use long.ld.missing = FALSE. The function will use SNPRelate::snpgdsLDpruning to prune the dataset. SNPs in LD are selected randomly. Default: long.ld.missing = FALSE.

  • ld.method: (optional, character) The values available are "composite", for LD composite measure, "r" for R coefficient (by EM algorithm assuming HWE, it could be negative), "r2" for r^2, "dprime" for D', "corr" for correlation coefficient. The method corr and composite are equivalent when SNPs are coded based on the presence of the alternate allele (0, 1, 2). Default: ld.method = "r2".

  • ld.figures: (logical) Generate long distance LD statistics and figures. Default: ld.figures = TRUE

  • path.folder: to write ouput in a specific path (used internally in radr). Default: path.folder = getwd(). If the supplied directory doesn't exist, it's created.

Interactive version

The workflow has a short-distance and an optional long-distance component:

  1. When loci contain multiple SNPs, choose the SNP-retention rule: 1 = mac, 2 = random, 3 = first, 4 = middle, or 5 = last.

  2. Answer "Do you want to continue filtering using long distance ld? (y/n):". Answering no returns after short-distance pruning.

  3. For reference-guided data, choose whether missingness should determine which SNP is retained from an LD group. De novo data currently uses basic SNPRelate pruning.

  4. Inspect the LD boxplot and answer "Enter the long LD threshold (filter.long.ld threshold, double/proportion):".

Use interactive.filter = FALSE and explicit short- and long-distance settings for a reproducible analysis. The missingness-aware option is passed as long.ld.missing through ....

References

Zheng X, Levine D, Shen J, Gogarten SM, Laurie C, Weir BS. (2012) A high-performance computing toolset for relatedness and principal component analysis of SNP data. Bioinformatics. 28: 3326-3328. doi:10.1093/bioinformatics/bts606

See also

SNPRelate

Author

Thierry Gosselin thierrygosselin@icloud.com

Examples

if (FALSE) { # \dontrun{
genome <- genometranslator::read_genome(
  data = "my.vcf",
  strata = "my.strata.tsv"
)

# Explore short- and long-distance LD interactively.
genome <- radr::filter_ld(data = genome)

# Alternatively, use a separate unfiltered GDS for scripted pruning.
scripted_genome <- genometranslator::read_genome("my_genome_scripted.gds")
scripted_genome <- radr::filter_ld(
  data = scripted_genome,
  interactive.filter = FALSE,
  filter.short.ld = "mac"
)

# Long-distance pruning can use missingness to choose among linked SNPs.
another_genome <- genometranslator::read_genome("my_genome_long_ld.gds")
another_genome <- radr::filter_ld(
  data = another_genome,
  interactive.filter = FALSE,
  filter.short.ld = "mac",
  filter.long.ld = 0.8,
  long.ld.missing = TRUE
)
} # }