Multi-Tissue eQTL code, part 1/4

Run Matrix eQTL for each tissue

Applying Matrix eQTL and recording the results for MT-eQTL analysis requies only minor adjustments to the standard code for running Matrix eQTL. Namely:

• Record all results by setting cis p-value threshold to 1:
     pvOutputThreshold.cis = 1;
• Turn off trans-eQTL discovery:
     pvOutputThreshold = 0;
• Turn off FDR calculation for memory efficiency and for faster processing.
     noFDRsaveMemory = TRUE;
• Record in a text file the number of degrees of freedom in the linear model for each tissue.
     cat(file="df.txt", tissue, "\t", me$param$dfFull, "\n", append=TRUE);

The code:

### Run separately for every tissue
tissue = "Adipose";

################ Running Matrix eQTL ##########################

library("MatrixEQTL");

### Load genotype info
snps = SlicedData$new();    
snps$LoadFile(paste0("genotypes/",tissue,".snps.txt"),
  skipRows = 1, skipColumns = 1, sliceSize = 500);

### Load gene expression info
expr = SlicedData$new();    
expr$LoadFile(paste0("expression/",tissue,".expr.txt"),
  skipRows = 1, skipColumns = 1, sliceSize = 500);

### Load covariates
cvrt = SlicedData$new();    
cvrt$LoadFile(paste0("covariates/",tissue,".covariates.txt"),
  skipRows = 1, skipColumns = 1, sliceSize = 500);

### Load gene locations
geneloc = read.table("geneloc.txt",
  sep = "\t", header = TRUE, stringsAsFactors = FALSE);

### Load SNP locations
snpsloc = read.table("snpsloc.txt",
  sep = "\t", header = TRUE, stringsAsFactors = FALSE);

options(MatrixEQTL.dont.preserve.gene.object = TRUE);

### Run Matrix eQTL
me = Matrix_eQTL_main(
  snps = snps,
  gene = expr,
  cvrt = cvrt,
  output_file_name = "",
  pvOutputThreshold = 0,
  useModel = modelLINEAR,
  errorCovariance = numeric(),
  verbose = TRUE,
  output_file_name.cis = paste0("eQTL_results_AL_",tissue,"_cis.txt"),
  pvOutputThreshold.cis = 1,
  snpspos = snpsloc,
  genepos = geneloc,
  cisDist = 1e5,      #### cis window = 100kb
  pvalue.hist = FALSE,
  noFDRsaveMemory = TRUE);

### Save the number of degrees of freedom for each tissue
cat(file="df.txt", tissue, "\t", me$param$dfFull, "\n", append=TRUE);


→ Next step: Combine results across tissues in one matrix of z-scores.

← Back to MT-eQTL page.

By Andrey Shabalin