Mercurial > repos > jasonxu > matrixeqtl

comparison MatrixEQTL/man/modelANOVA.Rd @ 0:cd4c8e4a4b5b draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Uploaded
author jasonxu
date Fri, 12 Mar 2021 08:12:46 +0000
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:cd4c8e4a4b5b
1 \name{modelANOVA}
2 \alias{modelANOVA}
3 \docType{data}
4 \title{
5 Constant for \code{\link{Matrix_eQTL_engine}}.
6 }
7 \description{
8 Set parameter \code{useModel = modelANOVA} in the call of \code{\link{Matrix_eQTL_main}} to indicate that the genotype should be treated as a categorical variable.
9 }
10 \note{
11 By default, the number of ANOVA categories is fixed to be 3.
12
13 To set it to a different number (say, 4) use the following command:
14
15 \code{options(MatrixEQTL.ANOVA.categories=4)}
16
17 To check the current settings run:
18
19 \code{getOption('MatrixEQTL.ANOVA.categories', 3);}
20 }
21 \examples{
22 library('MatrixEQTL')
23
24 # Number of columns (samples)
25 n = 100;
26
27 # Number of covariates
28 nc = 10;
29
30 # Generate the standard deviation of the noise
31 noise.std = 0.1 + rnorm(n)^2;
32
33 # Generate the covariates
34 cvrt.mat = 2 + matrix(rnorm(n*nc), ncol = nc);
35
36 # Generate the vectors with single genotype and expression variables
37 snps.mat = floor(runif(n, min = 0, max = 3));
38 gene.mat = 1 + (snps.mat==1) + cvrt.mat \%*\% rnorm(nc) + rnorm(n) * noise.std;
39
40 # Create 3 SlicedData objects for the analysis
41 snps1 = SlicedData$new( matrix( snps.mat, nrow = 1 ) );
42 gene1 = SlicedData$new( matrix( gene.mat, nrow = 1 ) );
43 cvrt1 = SlicedData$new( t(cvrt.mat) );
44
45 # name of temporary output file
46 filename = tempfile();
47
48 snps1
49 gene1
50
51 # Call the main analysis function
52 me = Matrix_eQTL_main(
53 snps = snps1,
54 gene = gene1,
55 cvrt = cvrt1,
56 output_file_name = filename,
57 pvOutputThreshold = 1,
58 useModel = modelANOVA,
59 errorCovariance = diag(noise.std^2),
60 verbose = TRUE,
61 pvalue.hist = FALSE );
62 # remove the output file
63 unlink( filename );
64
65 # Pull Matrix eQTL results - t-statistic and p-value
66
67 fstat = me$all$eqtls$statistic;
68 pvalue = me$all$eqtls$pvalue;
69 rez = c( Fstat = fstat, pvalue = pvalue)
70 # And compare to those from ANOVA in R
71 {
72 cat('\n\n Matrix eQTL: \n');
73 print(rez);
74 cat('\n R anova(lm()) output: \n')
75 lmodel = lm( gene.mat ~ cvrt.mat + factor(snps.mat), weights = 1/noise.std^2 );
76 lmout = anova( lmodel )[2, 4:5];
77 print( lmout )
78 }
79
80 # Results from Matrix eQTL and 'lm' must agree
81 stopifnot(all.equal(lmout, rez, check.attributes=FALSE))
82 }
83 \references{
84 The package website: \url{http://www.bios.unc.edu/research/genomic_software/Matrix_eQTL/}
85 }
86 \seealso{
87 See \code{\link{Matrix_eQTL_engine}} for reference and sample code.
88 }