In letaylor/variantBedOverlap: Calculates the overlap of variants in LD with features from BED files

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
# render to pdf
# R --slave -e 'rmarkdown::render("bioutils.Rmd", output_format = "rmarkdown::pdf_document")'
# render to html document specifically formatted for vignette
# R --slave -e 'rmarkdown::render("bioutils.Rmd", output_format = "rmarkdown::html_vignette")'

This document shows a demo of how to use variantBedOverlap. There is also a command line script.

```{bash cl, echo=T, eval=F, results='asis'}

get the lib dir for variantBedOverlap

install_dir=$(R --slave -e 'cat(find.package("variantBedOverlap"))')

see help options of command line script

Rscript "$install_dir/exec/variant_bed_overlap.R" --help

## (1) Get variants in LD

Get proxies from the 1000 Genomes Project via [proxysnps](https://github.com/slowkow/proxysnps).

```r
snps_q <- proxysnps::get_proxies(query = "rs2072014", pop = "FIN")
snps <- subset(snps_q, R.squared >= 0.8)

#knitr::kable( head(snps, 10) )
snps # enabled by setting --> df_print: paged

(2) Get BED overlaps

Get the overlaps of each variant (row) and genomic regions from a list of BED files. Here we load example BED files included in this package taken from Varsheny et al. 2017 (https://doi.org/10.1073/pnas.1621192114).

# few pre-packaged bed files from 
# https://theparkerlab.med.umich.edu/data/papers/doi/10.1073/pnas.1621192114/
dir <- system.file("extdata", package = "variantBedOverlap", mustWork = TRUE)

# get overlaps with all bed files in directory
snps_overlap <- variantBedOverlap::get_bed_overlaps(
  df = snps,
  dir = dir,
  col_itemRgb = 5
)

snps_overlap

(3) Plot the data

Plot the overlap data.

# xid_solid_line = list of x-axis IDs to add a line through
# varshney_chrhmm = flag to say assume BED file names are from Varshney et al
#                   2017. Given that assumption clean up the names to make them
#                   publication ready.
lst <- variantBedOverlap::plot_overlaps(
    df = snps_overlap,
    xid_solid_line = c("rs2072014", "rs35045598"),
    varshney_chrhmm = TRUE
)
print( lst$plt )

The output of plot_overlaps also contains the data underlying the plot. Note that ID is a factor now, sorted by POS. If varshney_chrhmm == TRUE, bed_feature will also be a factor sorted by chrhmm_state.

lst$df

letaylor/variantBedOverlap documentation built on Dec. 24, 2019, 11:21 a.m.