qtl2 package

knitr::opts_chunk$set(echo = TRUE)

:warning: this package is under development - use with care

`R/qtl2helper`

qtl2helper provides with some helper functions to work with objects created by several functions in the R/qtl2 package.

To install the package run:

# install the remotes package if not installed already
if (!("remotes" %in% rownames(installed.packages()))){
  install.packages("remotes")
}

# install the package directly from github
remotes::install_github("tavareshugo/qtl2helper")

Quick overview

Let's start with an example dataset from the qtl2 user guide:

library(qtl2)
library(qtl2helper)

# read example data
iron <- read_cross2( system.file("extdata", "iron.zip", package="qtl2") )
map <- insert_pseudomarkers(iron$gmap, step=1)
pr <- calc_genoprob(iron, map, error_prob=0.002)
Xcovar <- get_x_covar(iron)

# run QTL scan
out <- scan1(pr, iron$pheno, Xcovar=Xcovar)

# get QTL effects
c2eff <- scan1coef(pr[,"2"], iron$pheno[,"liver"])

qtl2helper provides methods for the tidy function (from the broom package), which convert several R/qtl2 objects to (long) data frames.

# convert scan1 object to a tibble
tidy(out)

You can provide a map of markers to tidy() to get marker positions in the output:

# convert scan1 object to a tibble
tidy(out, map)

This makes it ideal to do further plotting, for example with ggplot2:

library(ggplot2)

out_tbl <- tidy(out, map)

ggplot(out_tbl, aes(x = pos, y = LOD, colour = pheno)) +
  geom_line() +
  facet_grid(. ~ chrom, scales = "free_x", space = "free_x") +
  theme_classic() +
  theme(axis.text.x = element_blank(), axis.ticks.x = element_blank()) +
  labs(x = "Chromosome", y = "LOD score", colour = "Tissue")

Here's the same idea applied to the estimated QTL effects:

# tidy the scan1coef object
c2eff_tbl <- tidy(c2eff, map)

# plot
ggplot(c2eff_tbl, aes(x = pos, y = estimate, colour = coef)) +
  geom_line() +
  labs(x = "Chromosome 2 position", y = "QTL effects", colour = "Genotype")

Finally, here's an example of tidying permutations to obtain a tibble of thresholds:

# run permutation scans
operm <- scan1perm(pr, iron$pheno, Xcovar=Xcovar, n_perm=10)

# tidy output
operm_tbl <- tidy(operm)
operm_tbl

These can now be added to the scan plot, for example using geom_hline():

ggplot(out_tbl, aes(x = pos, y = LOD, colour = pheno)) +
  geom_line() +
  geom_hline(data = operm_tbl, 
             aes(yintercept = threshold, colour = pheno), 
             linetype = "dotted") +
  facet_grid(. ~ chrom, scales = "free_x", space = "free_x") +
  theme_classic() +
  theme(axis.text.x = element_blank(), axis.ticks.x = element_blank()) +
  labs(x = "Chromosome", y = "LOD score", colour = "Tissue")