knitr::opts_chunk$set( collapse = TRUE, comment = "#>" ) options(stringsAsFactors = FALSE)
episcan provides efficient methods to scan for pairwise epistasis in both case-control study and quantitative studies. It is suitable for genome-wide interaction studies (GWIS) by splitting the computation into manageable chunks. The epistasis methods used by episcan are adjusted from two published papers [@epiblaster; @epiHSIC].
# load package library(episcan)
First, we generate a small genotype dataset (
geno) with sample size of 100 subjects and 100 variables (e.g., SNPs) as well as a case-control phenotype (
set.seed(321) geno <- matrix(sample(0:2, size = 10000, replace = TRUE, prob = c(0.5, 0.3, 0.2)), ncol = 100) dimnames(geno) <- list(row = paste0("IND", 1:nrow(geno)), col = paste0("rs", 1:ncol(geno))) p <- c(rep(0, 60), rep(1, 40)) geno[1:5, 1:8]
To use episcan, simply start with the main function
episcan. There are three mandatory parameters to be set by the user: genotype data, phenotype data and phenotype category ("case-control" or "quantitative"). Since the data simulated above is not normalized yet, we need to set
scale = TRUE. By passing an integer number to parameter
chunksize, the genotype data will be split into several chunks of that size during the calculation. For the example above (using
chunksize = 20 means each chunk has 20 variables(variants) and the total number of chunks is 5. Moreover, in most cases, the result of epistasis analysis is huge due to the large number of the variable (variants) combinations. To reduce the size of the result file, setting a threshold of the statistical test (
zpthres) to have an output cut-off is a practical option.
episcan(geno1 = geno, pheno = p, phetype = "case-control", outfile = "episcan_1geno_cc", suffix = ".txt", zpthres = 0.9, chunksize = 20, scale = TRUE)
The result of
episcan is stored in the specified file ("episcan_1geno_cc.txt"). Let's take a look:
result <- read.table("episcan_1geno_cc.txt", header = TRUE, stringsAsFactors = FALSE) head(result)
In a genome-wide level epistasis study, it is usual to have millions of variables (variants). Analyzing such big data is super time-consuming. The common appoach is to parallelize the task and run the subtasks with High Performance Computing (HPC) techniques, e.g., on a cluster. By splitting genotype data per chromosome, the huge epistasis analysis task can be divided into relatively small tasks. If only 22 chromosomes exist in the initial task, there are 253 ((1+22)*22/2) subtasks after splitting and considering all the combinations of the chromosomes.
episcan supports two inputs of genotype data by simply passing information to "geno1" and "geno2".
# simulate data geno1 <- matrix(sample(0:2, size = 10000, replace = TRUE, prob = c(0.5, 0.3, 0.2)), ncol = 100) geno2 <- matrix(sample(0:2, size = 20000, replace = TRUE, prob = c(0.4, 0.3, 0.3)), ncol = 200) dimnames(geno1) <- list(row = paste0("IND", 1:nrow(geno1)), col = paste0("rs", 1:ncol(geno1))) dimnames(geno2) <- list(row = paste0("IND", 1:nrow(geno2)), col = paste0("exm", 1:ncol(geno2))) p <- rnorm(100) # scan epistasis episcan(geno1 = geno1, geno2 = geno2, pheno = p, phetype = "quantitative", outfile = "episcan_2geno_quant", suffix = ".txt", zpthres = 0.6, chunksize = 50, scale = TRUE)
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.