|
0
|
1 library("MatrixEQTL")
|
|
|
2
|
|
|
3 # Number of samples
|
|
|
4 n = 100;
|
|
|
5
|
|
|
6 # Number of variables
|
|
|
7 ngs = 2000;
|
|
|
8
|
|
|
9 # Common signal in all variables (population stratification)
|
|
|
10 pop = 0.2 * rnorm(n);
|
|
|
11
|
|
|
12 # data matrices
|
|
|
13 snps.mat = matrix(rnorm(n*ngs), ncol = ngs) + pop;
|
|
|
14 gene.mat = matrix(rnorm(n*ngs), ncol = ngs) + pop + snps.mat*((1:ngs)/ngs)^9/2;
|
|
|
15
|
|
|
16 # data objects for Matrix eQTL engine
|
|
|
17 snps1 = SlicedData$new( t( snps.mat ) );
|
|
|
18 gene1 = SlicedData$new( t( gene.mat ) );
|
|
|
19 cvrt1 = SlicedData$new();
|
|
|
20 rm(snps.mat, gene.mat);
|
|
|
21
|
|
|
22 # Slice data in blocks of 500 variables
|
|
|
23 snps1$ResliceCombined(500);
|
|
|
24 gene1$ResliceCombined(500);
|
|
|
25
|
|
|
26 # Produce no output files
|
|
|
27 filename = NULL; # tempfile()
|
|
|
28
|
|
|
29 # Perform analysis recording information for a histogram
|
|
|
30 me = Matrix_eQTL_main(
|
|
|
31 snps = snps1,
|
|
|
32 gene = gene1,
|
|
|
33 cvrt = cvrt1,
|
|
|
34 output_file_name = filename,
|
|
|
35 pvOutputThreshold = 1e-100,
|
|
|
36 useModel = modelLINEAR,
|
|
|
37 errorCovariance = numeric(),
|
|
|
38 verbose = TRUE,
|
|
|
39 pvalue.hist = 100);
|
|
|
40
|
|
|
41 # png(filename = "histogram.png", width = 650, height = 650);
|
|
|
42 plot(me, col="grey");
|
|
|
43 # dev.off();
|
